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numerical  equivalence  of  the  two  methods  have  been  demonstrated. 
It  should  be'noted  that  (whereas  in  the  potential  flow  problem 
viscosity  may  be  introduced  only  as  artificial  viscosity)  in  the 
rotational-flow  formulation  the  presence  of  viscosity  is 
consistent  with  the  formulation.  Therefore,  this  formulation  is 
applicable  to  the  solution  of  time-averaged  Na v i e r-Sto ke s 
equations:  in  particular,  a  simplified  thin-wake  analysis  with 

an  elementary  eddy-viscosity  model  for  turbulence  is  used  in  the 
numerical  applications.' 


Results  obtained  indicate  that,  for  appropriate  values  of 
emperical  parameters  (e.g.,  thickness  of  the  wake  at  the 
trailing  edge,  eddy  Reynolds  number,  intensity  of  the  far  wake 
sink  disk)  the  solution  reaches  steady-state,  with  the  wake 
converging  to  a  smooth  geometry  and  the  sectional  lift 
distribution  in  good  agreement  with  the  numerical  results  of  Rao 
and  Schatzle  and,  indirectly,  the  experimental  ones  of  Bartsch. 
Ground  effect  results  are  also  included.  However,  a  sensitivity 
analysis  over  the  effects  of  these  parameters  indicates  that  the 
results  were  not  always  as  good.  For  some  values  of  these 
parameters,  the  solution  may  reach  8  steady  state  but  the 
sectional  lift  distribution  does  not  necessarily  show  a  marked 
improvement  as  the  wake  roLls  up.  For  other  values  of  these 
parameters,  the  solution  may  be  completely  unstable.  The  main 
reason  for  this  sensitivity  to  the  values  of  the  empirical 
parameters  is  attributed  to  the  fact  that  concentrated  vortices 
are  used  for  the  discretization  of  the  formulation.  Additional 
work  is  needed  to  address  these  issues. 

^Analysis  of  the  computer  times  indicates  that  the  new 
formulation  presented  here  (for  rotational  flows)  has  much 
broader  applicability  than  the  formulation  for  potential  flows, 
while  requiring  approximately  the  same  amount  of  computer  time. 
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SECTION  1 


INTRODUCTION 


1.1  Motivation  and  Objectives 

Two  computational  methods  for  the  free  wake  aerodynamic  analysis  of 
a  helicopter  rotor  are  presented  in  this  report. 

The  availability  of  these  methods  considerably  enhances  the  computa¬ 
tional  capability  in  this  area.  This  capability  is  needed  in  many  aspects  of  the 
helicopter  design,  such  as  performance,  structural  analysis,  and  evaluation 
of  generalized  forces  for  flutter  analysis. 

In  all  of  these  cases  the  wake  geometry  is  an  essential  aspect  of  the 
problem.  The  approaches  used  in  the  past  can  be  classified  under  three 
broad  groups: 

A.  classical- wake  analysis 

B.  generalized-wake  analysis 

C.  free-wake  analysis 

In  the  classical  rotor-wake  formulation,  the  wake  for  a  hovering  rotor  is 
described  as  a  spiral  (helicoidal)  surface  which  is  obtained  from  the  assump¬ 
tion  of  uniform  vertical  flow.  This  method  is  not  only  limited  to  the  hovering 
case,  but  is  also  not  sufficiently  accurate  for  the  aerodynamic  analysis  of  he¬ 
licopter  rotors. 

On  the  other  hand,  a  generalized-wake  analysis  consists  in  an  analytical 
description  of  the  wake  geometry  in  terms  of  a  few  parameters  which  are 
evaluated  by  interpolation  of  experimental  results.  The  generalized-wake 
analysis  may  be  used  for  both  hover  and  forward  flight,  yields  accurate  re¬ 
sults  and  is  not  necessarily  more  expensive  than  the  classical-wake  analysis. 
However,  currently  it  requires  the  use  of  expensive  wind  tunnel  experiments 
for  the  generation  of  the  generalized-wake  model. 

Finally,  the  free-wake  analysis  consists  in  calculating,  step-by-step,  the 
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location  of  a  wake  point  at  a  given  time  step  from  its  location  at  the  preceding 
time  step:  the  drawback  with  this  approach  is  that  the  free-wake  analysis  is 
quite  expensive  in  terms  of  computer  time.  However,  it  could  be  used,  instead 
of  the  much  more  expensive  wind-tunnel  experiments  in  order  to  generate  a 
generalized  wake  model  for  hover,  forward-flight  and  maneuvering. 

1.2  Relation  of  Work  with  State  of  the  Art 

Two  methods  for  the  free-wake  analysis  of  rotors  in  incompressible  flows 
are  presented  here.  The  first  one,  for  potential  flows  is  an  extension  of  the 
work  by  Morino.  Kaprielian  and  Sipcic  (Refs.  1.  2  and  3).  The  wake  is 
treated  as  a  doublet  layer  (fully  equivalent  to  a  zero-thickness  vortex  layer). 
By  definition  of  potential  flows,  the  formulation  is  limited  to  zero-thickness 
wake.  [Artificial  viscosity  is  often  introduced  as  o.  numerical  approximation  to 
eliminate  numerical  instabilities,  but  is  difficult  to  justify  from  a  theoretical 
point  of  view  since  the  wake  is  treated  inconsistently:  as  zero-thickness  in  the 
-.  tegral  equation  for  the  evaluation  of  the  potential  on  the  surface,  ana  as 
finite-thickness  in  the  evaluation  of  the  velocity  field].  Therefore,  the  method 
is  not  expected  to  yield  good  results  (at  least  without  the  use  of  some  ad- 
hoc  numerical  scheme,  such  as  artificial  vorticity)  in  those  cases  for  which 
the  zero-thickness  assumption  is  too  drastic.  Examples  include:  a  rotor  in 
ground  effect  (which  is  the  specific  problem  addressed  here),  but  also  wake- 
fuselage  interaction,  blade-vortex  interaction,  blade  slapping  and  tip  vortex 
breakdown. 

The  second  method  presented  here  addresses  this  issue:  distributed  vor¬ 
ticity  is  allowed  in  the  flow  field.  In  this  case  the  flow  field  is  not  potential 
and  the  first  method  is  not  applicable.  The  second  method  is  based  on  the 
Helmholtz  decomposition  of  a  vector  field  into  a  solenoidal  part  and  al  ir- 
rotational  part.  Special  attention  is  given  in  this  report  to  the  relationship 
between  the  two  methods.  In  particular  it  is  shown  how  the  rotational- 
flow  approach  relates  to  the  potential- flow  approach  for  the  case  of  inviscid 
initially-irrotational  flows.  It  is  also  shown  how  the  method  can  be  extended 
to  study  viscous  flows:  this  yields  an  understanding  of  when,  why  and  how 
Euler  equations  are  a  good  approximation  to  Navier  Stokes  equations. 
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A  detailed  summary  of  this  report  is  given  in  Section  1.3  after  a  re¬ 
view  of  the  state  of  the  art.  We  feel  that  the  unsteady  analysis  is  the  main 
strength  of  the  methods  presented  here  over  other  existing  methods.  There¬ 
fore  particular  emphasis  is  given  in  this  review  to  the  unsteady-flow  analysis. 

An  excellent  review  on  aerodynamic  technology  for  advanced  rotorcraft 
was  presented  by  Landgrebe.  Moffitt.  and  Clark  in  Refs.  4  and  5.  Additional 
reviews  on  the  specific  issues  of  rotor  wakes  are  presented  in  Refs.  1  and  6-10. 
Therefore  only  works  which  are  particularly  relevant  to  the  objective  and  the 
motivation  of  this  work  are  included  in  this  brief  review,  which  is  not  to  be 
considered,  by  any  means,  complete.  Also,  only  works  dealing  with  free-wake 
analysis  of  rotors  in  incompressible  flows  are  considered  here.  (The  effects  of 
compressibility  are  dealt  with  in  a  separate  report  (Ref.  11)  and  therefore 
not  included  in  this  review'.  A  review  on  the  use  of  panel  methods  for  rotor 
aerodynamic  analysis,  also  not  included  here,  may  be  found  in  Ref.  1.) 

The  essence  of  the  state  of  the  art  in  incompressible  rotor  aerodynamics 
is  briefly  summarized  here.  As  mentioned  above,  the  various  aerodynamic 
analyses  of  the  rotor  fall  into  one  of  the  three  following  types: 

A.  Classical  wake,  i.e.,  a  wake  geometry  described  by  a  helicoidal  spiral 
with  the  pitch  obtained  by  assuming  uniform  flow. 

B.  Generalized  wake,  i.e.,  a  wake  geometry  obtained  by  interpolating 
experimental  data  in  terms  of  a  few'  parameters. 

C.  Free  wake.  i.e..  a  wake  geometry  obtained  computationally  as  an 
integral  part  of  the  solution. 

Generalized-wake  models  for  predicting  the  geometry  of  the  rotor  wake 
were  developed  from  experimental  data  by  Landgrebe  (Ref.  12  and  13). 
Crews,  Hohenemser  and  Ormiston  (Ref.  14)  and  Kocurek  (Ref.  15).  Land- 
grebe's  model  was  used  by  Rao  and  Schatzie  (Ref.  7)  in  their  lifting-surface 
theory:  their  work  is  particularly  significant  here  because  it  shows  that 
whereas  the  classical-wake  formulation  is  not  sufficiently  accurate,  a  con¬ 
siderable  improvement  in  the  comparison  with  experimental  results  of  Ref. 
16  can  be  obtained  simply  by  using  Landgrebe's  generalized-wake  geometry. 
Free-wake  analyses  are  considered  for  instance  by  Scully  (Ref.  17),  Pouradier 
and  Horowitz  (Ref.  18)  and  Summa  (Refs.  6.  19  and  20).  All  these  works  in¬ 
dicate  that  the  algorithms  used  are  unstable  unless  special  constraints  (such 
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as  specified  contraction  ratio  or  over-relaxation)  are  introduced.  Another 
important  issue  is  the  one  of  simplified  algorithms  which  can  be  used  for 
instance  to  model  the  far  wake:  several  models  are  available  for  the  hover 
case  (see  e.g..  Refs.  19-21).  However  none  of  these  models  is  applicable  to 
the  case  of  the  arbitrary  motion  considered  in  this  report.  The  fast-wake  for¬ 
mulation  which  Miller  (Ref.  21.)  introduced  for  modelling  the  whole  wake 
is  also  an  excellent  candidate  for  the  far-wake:  this  formulation  could  be 
extended  to  forward  flight  or  more  arbitrary  motion.  Additional  works  on 
w'ake  dynamics  which  should  be  mentioned  include  those  by  Kandil  (see  e.g.. 
Ref.  8)  and  Morino  and  Suciu  (Refs.  9  and  10).  The  formulation  of  Refs. 
9  and  10  is  extended  to  helicopter  rotors  in  Ref.  1  and  is  the  basis  of  the 
potential  methods  presented  here. 

1.3  Summary  of  Work 

As  mentioned  above,  this  report  includes  two  different  approaches  to 
study  the  flow  field  around  a  helicopter  rotor.  The  first  one  (see  Section  2) 
is  for  potential  flows  and  is  similar  to  that  of  Refs.  1-3.  [The  discussion  of 
the  boundary  condition  on  the  wake  and  of  the  trailing  edge  conditions  are 
more  rigorous  than  that  of  Refs.  1-3].  This  formulation  is  based  upon  an 
integral  equation  developed  by  Morino  (Refs.  22-240)  with  applications  to 
aircraft  (Refs.  25-27)  and  rotors  (Refs.  28-29).  The  wake  is  treated  as  a 
doublet  layer  (fully  equivalent  to  a  zero-thickness  vortex  layer). 

The  second  formulation,  presented  in  Section  3  is  for  incompressible  ro¬ 
tational  flows.  The  formulation  is  based  on  Helmholtz  scalar/vector-potential 
decomposition  (Refs.  30-32)  and  yields  the  following  computational  proce¬ 
dure:  given  the  vorticitv  distribution  at  time  t,  evaluate  the  velocity  ac¬ 
cording  to  Biot  and  Savart  law.  then  use  the  scalar-potential  contribution 
to  satisfy  the  normal  boundary  conditions,  and  finally  use  the  vorticitv  evo¬ 
lution  equation  to  evaluate  the  vorticitv  at  time  t  ~  dt.  A  novel  approach 
(based  on  the  use  of  materia!  curvilinear  coordinates)  is  used  to  solve  the 
vorticitv  evolution  equation  with  emphasis  on  a  thin  wake  approximation. 
The  formulation  is  presented  for  both  inviscid  and  viscous  flows  (with  a  very 
simple  eddy-viscositv  modelling  for  turbulence). 
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Section  4  deals  with  the  comparison  between  the  two  formulations.  In 
particular,  it  is  shown  how  in  the  limiting  case  of  zero-thickness  wake,  the 
rotational-flow  formulation  reduces  to  a  potential-flow  formulation  closely 
related  to  that  of  Section  2. 

Finally,  numerical  results  are  presented  in  Section  5.  and  conclusions  in 
Section  6. 


I 


I 
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SECTION  2 


POTENTIAL-FLOW  FORMULATION 

i 

A  general  formulation  for  the  analysis  of  potential  flows  is  presented  in 
this  section.  Special  emphasis  is  given  to  the  issues  related  to  the  wake 
transport,  which  were  first  introduced  in  Ref.  1.  This  section  includes  some 
refinements  with  respect  to  the  formulation  of  Ref.  1  (i.e..  the  formulation 
of  the  wake  condition  and  of  the  trailing  edge  conditions)  and  is  included 
here  essentially  for  the  sake  of  completeness.  The  differential  formulation 
(including  a  detailed  analysis  of  the  wake  condition)  is  presented  first  (Sec¬ 
tions  2.1  to  2.4).  followed  by  the  integral  formulation  (Sections  2.5  to  2.9). 

j  The  wake  generation  (trailing  edge  conditions)  is  discussed  in  Section  2.10, 

and  the  conditions  at  infinity  are  examined  in  Section  2.11. 

2.1  Governing  Equations 

|  In  this  report  we  will  assume  that  the  frame  of  reference  is  connected 

with  the  undisturbed  air.  We  also  assume  that  the  fluid  is  inviscid  and 
incompressible.  Hence,  the  phenomenon!  is  governed  by  the  continuity  equa¬ 
tion  for  incompressible  fluids  (conservation  of  mass) 


V  -v  =  0  (2  1) 

(where  v  is  the  velocity  vector)  and  by  Euler's  equations  (conservation  of 
momentum) 

~  =  -Vp  (2.2) 

Dt  p 

where  p  is  the  pressure,  p  the  density  (constant  for  incompressible  fluid),  and 
t  the  time,  whereas 

(25» 

is  the  material  or  substantial  derivative.  Equations  2.1  and  2.2  form  a  system 
of  four  partial  differential  equations  in  four  unknowns  vx ,  t>2,  us,  and  p. 
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In  order  to  complete  the  formulation  of  the  problem,  the  boundary  con¬ 
ditions  at  infinity,  on  the  body  and  on  the  wake  must  be  obtained.  Recalling 
that  the  frame  of  reference  has  been  assumed  to  be  connected  with  the  undis¬ 
turbed  air.  the  boundary  conditions  at  infinity  are 

v  =  o  (2.4) 


and 

p  =  Pc o  (2.5) 

On  the  surface  of  the  body  (rotor  in  our  case)  we  assume  that  the  surface 
is  impermeable.  This  implies  that  the  normal  components  of  the  velocity  v 
of  the  fluid,  and  of  the  velocity  vfc  of  the  body  coincide  at  point  x  on  the 
surface  ab  of  the  body: 


(v- vb) -11  =  0  (for  x  on  fffc)  (2.6) 

where  n  is  the  normal  to  ab  at  x. 

The  boundary  conditions  on  the  wake  are  discussed  in  Section  2.3  after 
introducing  the  concept  of  potential  wake. 

2.2  Wakes  in  Irrotational  Flows 

The  basis  of  the  discussion  on  wake  dynamics  is  the  well  known  Kelvin's 
theorem  which  states  that  the  circulation 

r  =  j,  v  •  dx  (2.7) 

over  a  material  contour  C  (i.e..  a  contour  which  is  made  up  of  material 
particles)  remains  constant  in  time  (<fr  dt  =  0).  As  shown  for  instance  in 
Ref.  30.  this  theorem  is  an  immediate  consequence  of  the  definition  of  T.  of 
Euler  equations  (Eq.  2.1)  and  of  the  fact  that  the  density  is  constant  (or,  in 
general,  that  the  fluid  is  barotropic). 

Next  assume  that  the  flow  field  is  irrotational  at  time  0.  Then  according 
to  Stokes  theorem 

(j}y  -dx  =  jj  V  >  v  •  n  d<r(x)  (2.8) 
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(where  C  is  the  contour  of  <r).  r  is  initially  equal  to  zero  for  any  contour  C 
connected  with  a  surface  a  fully  inside  the  fluid  volume.  Hence,  for  all  these 
contours,  r  remains  identically  equal  to  zero.  This  implies  that 

V  *  v  =  o  (2.9) 

for  almost  all  the  fluid  points  at  all  times:  the  only  points  to  be  excluded 
are  those  material  points  which  come  in  contact  with  the  solid  boundaries 
(since  a  contour  surrounding  any  of  these  points  in  general  is  not  completely 
inside  the  fluid  region  and  therefore  Kelvin’s  theorem  cannot  be  used  in  this 
case).  In  order  to  simplify  the  discussion  of  this  issue,  let  us  focus  on  the 
case  of  an  isolated  blade  with  a  sharp  trailing  edge  and  consider  only  those 
flows  for  which  the  fluid  leaves  the  surface  of  the  blade  at  the  trailing  edge 
(this  issue  is  further  discussed  later  in  Section  2.9).  These  flows  are  called 
attached  flows.  For  these  flows  the  points  which  come  in  contact  with  the 
rotor  are  only  those  emanating  from  the  trailing  edge  and  therefore  form  a 
surface:  such  a  surface  is  called  the  wake.  Equation  2.9  does  not  apply  to 
points  of  the  wake.  [Note  that  if  the  trailing  edge  is  fixed  with  respect  to 
the  frame  of  reference,  the  wake  is  a  streak  surface  (Ref.  32).  Since  the 
trailing  edge,  in  general,  moves  we  say  that  the  wake  is  a  generalized  streak 
surface.} 

We  may  conclude  that  for  an  inviscid  and  incompressible  fluid,  if  a  flow 
field  is  initially  irrotational,  it  remains  irrotational  at  all  times  except  for 
those  points  which  come  in  contact  with  the  surface  of  the  body.  If  the  flow 
is  attached  (i.e..  by  definition  if  the  fluid  leaves  the  body  only  at  the  points  of 
a  line  which  is  called  the  trailing  edge)  then  the  locus  of  these  points  forms 
a  generalized  streak  surface  which  is  called  a  wake. 

It  may  be  worth  emphasizing  that  the  result  of  zero-thickness  wake  is  a  di¬ 
rect  consequence  of  the  hypothesis  of  incompressible  (in  general,  barotropic) 
inviscid  flow  and  of  the  definition  of  attached  flows,  and  does  not  require 
additional  assumptions  (the  assumption  of  initially  irrotational  flows  is  not 
essential  for  the  thickness  of  the  wake  to  be  zero). 

2.3  Boundary  Conditions  on  Wake 

In  order  to  be  able  to  solve  the  mathematical  problem,  we  need  a  bound- 


ary  condition  on  the  wake.  This  condition  may  be  obtained  from  the  result 
shown  above  that  the  wake  has  zero  thickness,  and  hence,  is  a  surface  of  dis¬ 
continuity:  therefore,  the  principles  of  conservation  of  mass  and  momentum 
across  a  surface  of  discontinuity  may  be  used  to  obtain  the  wake  boundary 
conditions. 

Let  subscripts  u  and  l  indicate  the  two  sides  of  the  wake,  let  n  be  the 
outward  normal  on  side  u  and  let 

A/  =  /„  -  U  (210) 

denote  the  discontinuity  of  any  function  /  across  the  wake  surface.  For  the 
classical  wing- wake,  u  and  l  correspond  to  upper  and  lower  sides  respectively, 
n  is  the  upper  normal  and  A/  =  /upper  -  fiowtr-  In  the  case  of  a  rolled-up  rotor 
wake,  it  is  still  convenient  to  use  the  subscripts  u  and  i.  they  should  be 
understood  to  mean  the  sides  of  the  wake  that  are  the  continuation  of  the 
upper  and  lower  sides  of  the  blade  respectively.] 

The  equations  of  conservation  of  mass  and  momentum  across  a  possible 
surface  of  discontinuity  (e.g..  a  shock  wave  or  a  wake)  are  given  by  (Ref.  30) 

A  p(t> n  -  V.)  -  0  (2.11) 

A  p{  vn  -  t,)v  -  pn  =  />(  tn  -  i-,)Av  *  Ap  n  =  0  (2.12) 

where  vn  =  v  •  n  is  the  normal  component  of  the  velocity  v,  and  v,  is  the 
velocity  of  the  surface  (by  definition,  in  the  direction  of  the  normal  ■).  Note 
that 

Av,  =  0  (2.13) 

In  the  case  of  incompressible  flows,  Ap  =  0  and  hence,  Eqs.  2.11  and  2.13 
yield 

At'„  =  0  (2.14) 

and  the  only  possible  surfaces  of  discontinuity  are  waves.  (Surfaces  of  discon¬ 
tinuity  with  Ai„  »  0  are  called  shocks,  and  are  not  possible  in  incompressible 
flows).  Then  the  normal  component  of  Eq.  2.12  yields 


Hence  using  Eq.  2.15,  Eq.  2.12  reduces  to 


p(vn  -  u,)Av  =  o  (2.16) 

This  implies  either  Av  =  0  (i.e.,  the  surface  under  consideration  is  not  a 
surface  of  discontinuity)  or.  if  there  is  a  discontinuity  (i.e.,  on  a  wake  surface), 

un  =  t,  (2.17) 

which  indicates  that  the  fluid  does  not  penetrate  the  surface  of  the  wake. 


2.4  Potential  Flow;  Potential  Wakes 


Next  consider  a  well  known  theorem  from  vector  field  theory  which  states 
that  if  a  vector  field  v  is  irrotational.  then  there  exists  a  function.  p.  called 
velocity  potential  such  that 

V  =  tv  (2.18) 

If  the  region  is  simply  connected  (as  in  the  case  considered  here),  then  p  is 
single- valued. 

Hence,  the  results  of  Section  2.3  may  be  restated  as  follows:  for  an 
inviscid  incompressible  fluid,  a  flow  field  which  is  attached  and  initially  ir¬ 
rotational  is  potential  at  all  times  and  at  all  points,  except  at  the  points  of 
wake. 

If  the  flow  is  potential,  i.e..  if  v  is  given  by  Eq.  2.18.  Eq.  2.1  may  be 
rewritten  as 

V2p  =  0  (2.19) 


where  V2  is  the  Laplacian  operator. 

Similarly.  Eq.  2.4  may  be  rewritten  as 


p  —  0  (for  x  ai  sc)  ( 2 . 20 ) 

whereas  Eq.  2.6  becomes 

dp 

—  =  v>  a  (for  x  on  crb )  (2.21) 

an 

Furthermore.  Eq.  2.2  may  be  integrated  to  yield  Bernoulli’s  theorem 


(taking  the  gradient  of  Eq.  2.22  and  using  Eq.  2.18  one  obtains  Eq.  2.2:  the 
constant  on  the  right  hand  side  of  Eq.  2.22  is  obtained  from  the  boundary 
conditions  at  infinity.) 

Next,  in  order  to  obtain  a  convenient  condition  on  the  wake,  combine 
Bernoulli's  theorem  (Eq.  2.22)  with  the  wake  condition.  Ap  =  0  (Eq.  2.15). 
This  yields 


dt 


dJEL  _  I(v 

dt  V  “ 


v,  •  V, )  =  0 


(2.23) 


that  is. 


dpu 

dt 


dpi 

dt 


.  .  ,  d<Pu  dpi 

■  V( )  ■  Vu  -  V/ )  =  -  -  - - 

1  '  '  dt  dt 


(vu  -  v,)  •  Av  =  0 


(2.24) 


[Note  that  Eq.  2.14  implies  in  particular  that  only  the  tangential  com¬ 
ponents  of  the  velocity  are  discontinuous  across  the  surface  of  the  wake,  that 
is  (using  local  coordinates  ,  p  and  £s  with  p  and  p  on  the  surface  of  the 
wake  and  p  in  the  direction  of  the  normal  n) 


Av  =  Aft’jg1  -  v2g2  -  vsn)  =  Aft^g1  +  At>2g2) 

=  A  (  **  g‘  -  %2)  =  A^g>  -  A^g 
l 8  dp*  j  dp*  dp* 


d? 
d&<p 


d^ip  j 

=  l^g  *  ~dp 


dp 

g3=  V,(A*>) 


(2.25) 


where  va  are  the  covariant  components  of  v.  whereas  g“  are  the  contravariant 
base  vector,  related  to  the  covariant  base  vector 
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by  the  relationship  ga  g3  =  6%  where  6*  is  the  Kronecker  delta.  In  addition 

Vtf  =  IfT*1  *  ^«2  <2-27> 

is  the  tangential  portion  of  the  vector  V  f.  It  is  important  to  note  that  V(A <p) 
is  a  meaningless  expression  since  A p  is  defined  only  over  the  surface  of  the 
wake.] 

Also,  it  is  possible  to  write 


dpu  _  d<p,  _  /  d  d  \ 

dt  dt  \dt~Vndnl  * 


(2.28) 
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The  right  hand  side  of  Eq.  2.28  is  meaningful  (time  derivative  for  an  observer 
moving  with  the  surface  in  the  direction  of  the  normal  to  the  surface) ,  even 
though  d&p/dt  and  d&p/dn  individually  are  not. 

Therefore,  introducing  the  concept  of  a  wake  point  xw  as  a  point  having 
velocity 


(2.29) 


and  using  Eqs.  2.25  and  2.28.  Eq.  2.24  may  be  rewritten  as 


Dw 

~Dt 


M  =  o 


(2.30) 


(2.31) 


where 

Du,  _  d  _  _ 

Dt  ~  dt  Vw‘  dt*^dC  ~wd?  ~n  dn 
is  the  material  time  derivative  for  a  function  defined  only  over  a  surface. 
Equation  2.2  may  be  integrated  to  yield 


\p  =  constant  in  time  (2.32) 

following  a  wake  point  xu  which  has  velocity  vw  given  by  Equation  2.29.  The 
value  of  Ap  is  equal  to  the  value  that  it  had  when  the  point  xw  left  the  trailing 
edge. 

Note  that  t>*  and  vj  are  the  average  values  of  the  contravariant  tangential 
components  of  v„  and  v(  (recall  that  the  tangential  components  of  v  are 
discontinuous  across  the  wake),  whereas  u*  =  vn  =  v,  (see  Eq.  2.17)  since 
vn  is  continuous  across  the  wake.  It  may  be  worth  noting  that  nj  does  not 
appear  in  Equation  2.24.  However  the  fact  that  \w  has  normal  velocity  equal 
to  vn  stems  from  the  use  of  Equation  2.28. 


2.5  Summary  of  Differential  Formulation 

In  summary,  if  the  flow  field  of  an  inviscid  incompressible  fluid  is  initially 
irrotational.  it  remains  irrotational  at  all  times  except  for  the  points  of  the 
wake,  and  therefore  except  for  those  points,  may  be  expressed  as 

v  =  Vp  (2.33) 

Then  (see  Eq.  2.11),  p  satisfies  Laplace's  equation 

VV  =  0  (2.34) 
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with  condition  at  infinity  (see  Eq.  2.20) 


-  0 


(2.35) 


whereas  the  boundary  condition  on  the  rotor  blade  is  given  by  (see  Eq.  2.21) 


^  =  vt  n  (2.36) 

On  the  other  hand,  for  points  on  the  wake,  the  velocity  is  given  by  (see  Eq. 
2.29) 

v=^(V^u-V^)  (2.37) 

whereas,  (see  Eq.  2.32) 

=  constant  in  time  (2.38) 


following  x„, . 

Equations  2.34  to  2.38.  with  the  addition  of  the  trailing-edge  condition 
(see  Section  2.10),  may  be  used  to  obtain  the  solution  for  p.  Once  p  is  known, 
the  velocity  may  be  evaluated  using  Eq.  2.33.  Then  the  pressure  may  be 
evaluated  using  Bernoulli's  theorem  (see  Eq.  2.22) 


___  1  ,2  P  _  Poo 

dt  2  p  p 


(2.39) 


2.6  Green's  Theorem  for  Laplacian  Operator 

The  integral  equation  used  in  this  work  is  a  particular  case  of  that  in¬ 
troduced  by  Morino  (Refs.  22-24:  see  also  Refs.  25-29  for  applications  to 
airplanes  and  rotors)  for  the  general  case  of  potential  compressible  flows 
for  bodies  having  arbitrary  shapes  and  motions.  The  derivation  of  such  an 
equation  is  presented  here  for  the  specific  case  of  interest,  incompressible 
potential  flow  around  a  helicopter  rotor  having  arbitrary  shape  and  motion. 
The  integral  equation  is  based  on  the  Green's  function  method,  outlined  in 
the  following. 

The  basis  for  the  method  used  here  is  the  second  Green's  identity 

n  ,0) 


where  /  =  /( y).  g  =  g(y).  and 


Equation  2.40  is  easily  obtained  by  combining  Gauss  theorem. 

jJJ  a  iV(y)  =  -  Jf  a •  n  da(y) 


(2.41) 


(2  42) 


where  a  surrounds  V  and  n  is  the  inward  directed  normal,  with  the  identity 


*vU*yg-g*vf)  =  f*lg-gV2vf  (2-43) 

The  essence  of  the  method  consists  in  choosing  for  g  the  Green's  function, 
G.  defined  by 

V*G  =  *(y-x)  (2.44) 

where  6  is  the  Dirac  delta  function  defined  by 

jjpJ(y)S(y)dV(y)  =  f{Q)  (2  45) 

The  boundary  condition  at  infinity  for  Eq.  2.44  is 


G  =  0 


(2.46) 


The  solution  to  Eqs.  2.44  and  2.46  is  the  well  known  potential  for  the  unit 
source 


G  =  - 


4xr 


where 


r  =  y  -  * 


(2.47) 

(2.48) 


Introducing  the  function 


I 


a 


i 

i 


E( y)  =1  y  inside  a 

=  0  y  outside  a 


Equation  2.46  implies 

f£  f\y)Hy)dV(y)  =  jjj°°  E(y)f(y)6(y)dV(y)  =  £7(0)/( 0) 


(2.49) 


(2.50) 
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Combining  Eqs.  2.40.  2.44.  2.47  and  2.50  yields  the  elegant  formula 


£W/W  *  £  \jTn  TTr  ~  4  (is)j  d°M- SSL  ^  ^  dVM 


(2.51) 


Equation  2.51.  the  key  to  the  method  presented  here,  does  not  usually  have 
a  name  associated  with  it  (except  as  being  the  key  to  Greens 's  function 
method):  we  will  refer  to  it  as  Green's  theorem  for  the  Laplacian  operator. 


2.7  Green's  Theorem  for  Incompressible  Potential  Aerodynam¬ 


ics 


In  order  to  obtain  Green's  theorem  for  incompressible  aerodynamics,  we 
identify  the  function  /  in  Eq.  2.51  with  the  velocity  potential.  Hence  the 
volume  V  must  exclude  the  volume  V6  of  the  rotor  as  well  as  a  thin  layer, 
Vw.  which  includes  the  wake  surface  uw  (because  the  equation  V2y?  =  o  is  not 
valid  in  Vb  and  on  aw).  Finally,  the  volume  V  must  be  bound  by  an  outer 
surface  Hence  the  surface  a  (the  boundary  of  the  volume  V )  is  formed 
by  two  surfaces:  the  surface  <rbu,  (which  surrounds  the  rotor  volume,  Vb,  and 
the  wake  volume,  Vw:  note  that  the  inward  normal  for  a  is  outward  for  <rkw), 
and  the  surface  <rx  (which  is  a  spherical  surface  of  radius  R  and  center  x). 

Then,  using  Eq.  2.34.  Eq.  2.51.  for  /  =  p,  reduces  to 


da(y) 


(2.52) 


Note  that,  under  the  conditions  [it  is  verified  a  posteriori  (Eqs.  2.79  -  2.80) 
that  these  conditions  are  satisfied] 


lim  =  0 
a— o© 


df> 

lim  /?  — —  =  0 


a— oc  dn 

as  the  radius  R  of  goes  to  infinity,  one  obtains 

lim  ff  fi-,  i(fi)l,to-0 

a  —  oc  JJ 0 [cm  4irr  dn  V4rr/j 

Therefore,  as  R  tends  to  infinity.  Eq.  2.52  becomes 


(2.53) 

(2.54) 


(2.55) 


dy  - 1 

d 

( zi  V 

dn  4xr 

r  dn 

V  4  Trr  / 

d<r(y) 


(2.56) 
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\  V  *  %  % 


V  k,\  \ 


Next  let  the  two  sides  of  the  surface  <r'a  surrounding  the  wake  become 
infinitesimally  close  to  the  surface  of  the  wake.  Note  that  in  this  process,  the 
closed  surface  o'u  of  the  wake  is  replaced  by  the  two  sides  of  an  open  surface 
<tu,  .  Let  n  be  the  normal  on  the  side  u  of  aw.  In  the  limit,  one  obtains 


JJ0'w  dn  4 vr  \  dn  j  4 vr 

(since,  see  Eq.  2.14.  A(dyr/d«)  =  Avn  =  0).  whereas 


(2.57 


d_ 

dn 


(2.58) 


Therefore,  in  the  limit.  Eq.  2.52  reduces  to 
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where  ah  is  the  (closed)  surface  of  the  rotor  blade  and  aw  is  the  (open)  surface 
of  the  wake  of  the  rotor  blade.  Furthermore  (see  Eq.  2.14),  Ayr  =  ysu  -  v>t, 
where  n  is  the  normal  on  the  side  u  of  the  wake.  Note  that  A <p  on  cw  is 
evaluated  from  Eq.  2.38.  Note  also  that  the  vortex-layer  wake  of  the  rotor 
is  represented  as  a  doublet  layer.  The  proof  of  the  equivalence  of  doublet 
layers  and  vortex  layers  is  given,  for  instance,  in  Reference  31. 


2.8  Integral  Equation 


Equation  2.54  may  be  used  to  obtain  the  value  of  #  at  any  point  in  the 
field  if  the  value  of  ys  and  d<pjdn  on  the  surface  of  the  rotor  and  the  value 
of  Ayr  on  its  wake  are  known.  Note  that  d<pu dn  is  known  from  the  boundary 
condition  on  the  rotor,  Eq.  2.36.  whereas  Ayr  may  be  calculated  from  the 
boundary  condition  on  the  surface  of  the  wake.  Eq.  2.38. 

Hence,  in  order  to  be  able  to  use  Eq.  2.59.  one  must  have  an  equation 
for  evaluating  p  on  the  surface:  such  an  equation  is  obtained  by  noting  that 
if  x  approaches  a  point  on  the  surface  of  the  rotor,  then  the  value  of  »r(x) 
approaches  the  value  of  ^  on  the  surface  of  the  rotor.  In  order  to  perform 
this  limit,  it  is  convenient  to  interpret  the  doublet  integral  in  terms  of  solid 
angles.  Note  that  (see  Figure  1). 


d 

dn 


n 


da{y) 
F  4  z-rS 


1 


4r 


cos  a 


<My) 

r2 


^rfn(y) 

4  n 


(2.60) 


2  -  11 


where  r  =  y  -  x,  and 


jnr  >  (My)  cos  a 

rffl(y)  =  — ^ — 

is  the  elementary  solid  angle,  i.e..  the  surface  element  that  da  projects  upon 
the  unit  sphere.  Hence 

Note  that,  as  is  well  known 

fi(x)  =  4ir  for  x  inside  <rb 

=  2ir  for  x  on  ab  (regular  point) 

=  0  for  x  outside  <7s  (2.62) 

(where  a  regular  point  is  a  point  where  ab  has  a  unique  tangent  plane). 
Comparing  Eqs.  2.46  and  2.62  yields,  for  x  both  inside  and  outside  a. 

£(x)  =  t  -  ^  (2.63) 

Using  Eqs.  2.61  and  2.63,  Eq.  2.53  may  be  rewritten  as 

The  advantage  of  Equation  2.64  over  Equation  2.59  is  that  each  term  is 
continuous  (as  x  crosses  ab)  and  therefore  is  valid  in  particular  even  if  x  is  on 
ab.  This  implies  that  Eq.  2.63  is  also  valid  even  if  x  is  on  ab  (whether  x  is  a 
regular  point  of  ab  or  not).  In  particular,  if  x  is  a  regular  point  of  ab.  then 
n  =  2*-  and  E  =  1/2.  Equation  2.47  may  thus  be  generalized  as 

E(x)  =  1  -  H(x)  4;r  —  1  x  outside  ab 

=  -  x  on  ab  (regular  point) 

-  0  x  inside  ab  (2.65) 


In  any  event,  for  x  on  ah.  Eq.  2.59  (with  E  given  by  Eq.  2.65)  is  an  integral 
equation  relating  the  unknown  values  of  the  velocity  potential  on  the  surface 


of  the  rotor,  to  the  values  of  dp'dn  (prescribed  by  the  boundary  condition 
on  the  surface  of  the  blade)  and  the  values  of  the  potential  discontinuity  on 
the  wake  (known  from  the  preceding  time  history). 

2.9  Velocity  Field:  Wake  Transport 


Once  p  on  the  surface  is  known.  Eq.  2.59  may  be  used  to  calculate  p 
anywhere  in  the  field,  and  hence  the  velocity  at  any  point  in  the  field  as.  (see 
Eq.  2.33) 


V  = 

=  ** 

(2.66) 

Noting  that 

'0 

K 

(2  67) 

where 

r  = 

y  -  x 

(2  68) 

one  obtains 

dp  -r  d  - 

dn  4 ffr3  ^  dn 

(2.69) 

This  equation  may  also  be  used  to  calculate  the  velocity  of  the  points  of  the 
wake  if  the  contribution  of  the  wake  integral  in  an  infinitesimal  neighborhood 
of  the  wake  point  is  excluded  from  the  calculation:  such  a  contribution  is 
responsible  for  the  velocity  discontinuity  across  the  wake  and  its  exclusion 
automatically  yields  the  semi-sum  between  the  two  values  on  the  two  sides 
of  the  wake. 

2.10  Wake  Generation:  Trailing  Edge  Condition 

In  this  section  the  issue  of  the  wake  generation  is  analyzed.  First  the 
analytical  results  of  Section  2.8  are  interpreted  from  a  physical  point  of  view. 
Next,  the  trailing  edge  condition  is  introduced,  and  then  the  mechanism  of 
wake  generation  is  discussed. 

In  order  to  discuss  the  problem  of  wake  generation,  it  is  convenient  to 
consider  the  problem  of  a  rotor  subject  to  a  sudden  start,  i.e..  the  case  of  a 
rotor  which  for  t  <-  o  is  at  rest,  and  is  surrounded  by  a  fluid  which  is  also  a 
rest.  At  time  o~.  the  rotor  is  assumed  to  have  finite  velocity. 
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Since  for  time  t  <  0.  the  fluid  was  at  rest,  the  wake  has  not  yet  had  time 
to  develop.  Therefore  A*;  =  0  everywhere  in  the  field  and  the  solution  to  the 
problem  at  time  0~  is  obtained  by  solving  Equation  2.59  with  Ap  =  0.  that  is 


£(x)*r(x) 


dp  1 

d 

f  JLV 

dn  4v t 

i  4  ffr  / 

(2.70) 


with  E  given  by  Equation  2.65. 

l  This  corresponds  to  the  solution  of  the  Laplace  equation  without  wake 

contribution  and,  as  well  known  (for  instance  from  conformal  mapping  solu¬ 
tion  in  two  dimensional  flows),  this  yields  a  separation  point  which  in  general 
is  different  from  the  trailing  edge,  and  hence  causes  infinite  velocity  at  the 
trailing  edge.  In  turn,  this  indicates  the  presence  of  a  vortex  at  the  trailing 
edge:  this  manifests  itself  by  the  fact  that  Apu  ^  0. 

It  is  known  that  (since  our  potential  is  single  valued)  the  solution  to  the 
exterior  Neumann  problem  for  Laplace  equation  is  unique.  (Multi-valued 
potential  functions  occur  in  multiply-connected  regions). 

[The  proof  is  easily  obtained,  (see,  e.g.  Ref.  32).  by  integrating  *\e 
identity 

V  •  (uVu)  =:  uV2u  —  Vu  ■  Va  (2.71) 

to  obtain 

=  12721 

The  difference,  u,  between  two  solutions  to  the  Neumann  problem  satisfies  the 
Laplace  equation  with  homogeneous  boundary  conditions  on  a.  Therefore, 
the  first  two  integrals  in  Eq.  2.72  are  equal  to  zero.  This  implies  Vu  =  0.  i.e.. 
the  difference,  u.  between  two  solution  is  a  constant,  which  is  equal  to  zero 
because  of  the  condition  at  infinity:  in  other  words  the  solution  is  unique.] 
This  is  in  sharp  conflict  with  the  classical  two-dimensional  result  in  which 
a  suitable  vortex  can  be  added  inside  an  airfoil  to  obtain  smooth  flow  at 
the  trailing  edge.  The  mathematical  reason  for  the  difference  is  that  the 
flow  region  around  an  airfoil  is  doubly  connected,  whereas  that  around  a 
rotor  (as  well  as  most  three  dimensional  shapes  of  aerodynamic  .interest)  is 
simply  connected.  From  an  intuitive  point  of  view,  we  can  say  that  in  a 
doubly  connected  region,  it  is  possible  to  add  a  vortex  of  arbitrary  intensity 


inside  an  airfoil  (or  inside  a  doughnut-shaped  object,  which  also  corresponds 
to  a  doubly-connected  flow-region).  Then,  adding  a  suitable  single- valued 
solution,  one  may  obtain  a  nontrivial  solution  to  the  homogeneous  Neumann 
exterior  problem  for  the  Laplace  equation.  This  solution  can  always  be  added 
to  the  solution  of  the  airfoil  problem,  which  is  therefore,  not  unique.  It  should 
be  noted  that  the  potential  corresponding  to  a  vortex  is  multivalued  and  the 
above  proof  of  uniqueness  fails  for  multivalued  potential  functions. 

In  the  case  of  a  rotor,  it  is  possible  to  add  a  vortex  inside  the  rotor- 
surface:  however,  it  is  not  possible  to  have  a  contour  which  is  ‘interlocked' 
with  the  vortex  without  penetrating  the  rotor  surface  (this  is  true,  by  defi¬ 
nition,  for  any  simply  connected  region).  The  potential  flow  is  then  single¬ 
valued  and  is  impossible  to  generate  a  nontrivial  solution.  The  solution  to 
the  rotor  problem  is  thus  unique. 

In  order  to  clarify  the  issue  of  the  trailing  edge  conditions,  we  will  call 
the  Kutta-Joukowski  condition  as  that  smooth  flow  trailing  edee  condition 
which  is  used  to  eliminate  the  issue  of  a  nonunique  solution  (see  Refs.  33 
and  34).  This  condition  is  required  for  two-dimensional  steady  flows,  but 
not  in  three-dimensional  flows,  unless  the  region  is  multiply  connected,  nor, 
for  that  matter,  for  two  dimensional  unsteady  flows. 

In  unsteady  flows,  the  wake  generation  is  responsible  for  the  elimination 
of  singularities  at  the  trailing  edge.  The  correct  trailing  edge  condition  for 
unsteady  flows,  we  believe,  that  concentrated  vortices  typically  do  not  form 
at  the  trailing  edge. 

More  precisely,  we  introduce  the  following  assumption:  if,  at  certain 
discrete  times,  i  =  T, ,  concentrated  vortices  form  at  the  trailing  edge  (this  may 
happen  because  of  any  time-discontinuity  in  the  velocity  distribution  on  the 
surface  of  the  rotor,  i.e.,  infinite  acceleration  as  in  the  case  of  the  sudden  start 
discussed  above),  then  these  vortices  are  immediately  shed.  (This  condition 
must  be  satisfied  if  we  want  the  solution  to  our  problem  to  be  the  limiting 
case  of  the  solution  of  Navier-Stokes  equations).  At  any  other  time  (when  the 
acceleration  is  finite),  there  are  no  concentrated  vortices  at  the  trailing  edge. 
The  implication  of  this  assumption  is  that  (except  at.f  =  r,)  the  value  of  Ay? 
is  continuous  at  the  trailing  edge:  for  the  well  known  equivalence  between 
doublet  and  vortex  layers  (Ref.  31)  indicates  a  discontinuity  in  the  doublet 
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intensity  implies  the  presence  of  a  concentrated  vortex.  In  other  words,  if  xu 
is  a  point  on  the  wake  and  x,,  is  a  point  on  the  trailing  edge. 


]im  A^>(Xu,)  =  A  fit  (2.73) 

where  A^t  is  the  difference  in  the  values  of  the  potential  as  the  points  on 
upper  and  lower  surfaces  of  the  rotor  blade  approach  a  point  of  the  trailing 
edge. 

At  the  times  t  =  T,  where  the  acceleration  is  infinite.  Eq.  2.73  may  be 


generalized  as 

lim  Av?u,(t^  -r  «) 

-  AvJ,e(C) 

(2.74) 

and 

lim  (t~  -  «) 

<  — »0 

=  A<p,e(t") 

(2.75) 

which  expresses  that  the  intensity  of  A*>  on  the  wake  has  a  discontinuity. 
The  implication  is  that  a  vortex  is  shed  at  t  =  T,.  and  is  convected  with  the 
wake. 

Note  that  without  this  asumption.  the  solution  si  not  unique.;  there  exist 
at  least  two  solutions:  that  considered  here,  and  another  in  which  the  vortex 
remains  at  the  trailing  edge.  It  would  be  desirable  to  prove  that  under  the 
trailing  edge  condition  proposed,  there  exists  a  theorem  of  existence  and 
uniqueness  for  the  equation  of  potential  flows. 


2.11  Conditions  at  Infinity 


In  this  subsection.  Eqs.  2.53  and  2.54  are  verified.  Note  as  r  goes  to 
infinity.  Equation  2.56  yields 


(2.76) 


Even  if  the  flux  through  <rbw  is  not  equal  to  zero.  i.e.. 

$  * 0 


(2  77) 
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SECTION  3 


ROTATIONAL-FLOW  FORMULATION 


The  formulation  for  incompressible  rotational  flows  is  presented  in  this 
Section.  Following  the  classical  approach  of  Prandtl.  the  flow  region  is  di¬ 
vided  into  two  regions:  an  inner  one  (i.e..  the  boundary  layer  region,  near 
the  boundary  surface)  and  an  outer  one.  The  emphasis  is  on  the  outer  region 
(the  inner  one  is  assumed  to  be  solved  by  a  boundary  layer  code  or.  if  neces¬ 
sary.  by  a  thin-layer  Navier-Stokes  formulation).  The  formulation  is  based  on 
Helmholtz  scalar/vector-potential  decomposition.  Whereas  the  derivation  is 
strictly  mathematical,  the  resulting  formulation  may  be  restated  in  physical 
terms  as  follows.  If  the  vorticity  at  time  t  is  known,  then,  using  Biot  and 
Savart  law  (Section  3.3),  one  obtains  the  velocity  •‘induced  by  the  vorticity”; 
then  the  scalar-potential  contribution  is  added  to  satisfy  the  normal  bound¬ 
ary  condition  between  inner  and  outer  region.  The  two  above  steps  allow  to 
evaluate  the  velocity  field  from  the  vorticity  field.  Once  the  velocity  field  at 
time  t  is  known,  the  loop  may  be  closed  by  evaluating  the  vorticity  field  at 
time  t  -  dt  (from  that  at  time  t)  by  using  the  vorticity  evolution  equation. 

A  new  approach  for  solving  the  vorticity-evolution  equation  is  used  here. 
For  the  sake  of  clarity,  the  formulation  for  inviscid  flows  is  presented  first  (the 
formulation  is  extended  to  include  diffusion,  via  a  very  simple  eddy-viscosity 
turbulence  modelling,  in  Section  3.9).  As  shown  by  Morino  (Refs.  35  and 
36),  for  inviscid  flows,  the  material  contravariant  components  (i.e..  the  con- 
travariant  components  in  a  curvilinear  coordinate  system  wrhich  moves  with 
the  material  points.  Section  3.5)  of  the  vorticity  remain  constant  in  time 
(Section  3.6).  In  particular,  for  thin-wake  approximation,  this  implies  that 
the  wake-thickness  remains  constant  (Section  3.7).  This  result  is  used  for 
the  numerical  solution  of  the  vorticity  evolution  equation:  given  the  velocity 
field,  one  may  evaluate  the  new  location  of  the  wake:  from  this,  one  evaluates 
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the  base  vectors  and  the  new  values  of  the  vorticitv.  as  a  linear  combination 
of  the  (time-independent)  material  contravariant  components  of  the  vorticity 
and  their  corresponding  base  vectors  (see  Eq.  3.35).  As  mentioned  above,  the 
formulation  is  then  extended  to  the  ensemble-averaged  Navier-Stokes  equa¬ 
tions  (with  a  simple  eddy-viscosity  modelling  for  turbulence):  the  evolution 
equation  for  the  contravariant  material  components  is  obtained  and  used  to 
derive  the  expression  for  the  time-evolution  of  the  wake  thickness.  [Actual 
numerical  calculations  indicate  that  the  effect  of  the  thickness  change  may 
be  important  in  avoiding  numerical  instabilities  (see  Section  5)]. 

3.1  Governing  Equations 

The  assumption  of  incompressible  flow  is  used  throughout  this  report. 
However,  it  is  important  to  differentiate  between  incompressible  fluids  and 
incompressible  flows.  An  incompressible  fluid  has  zero  compressibility,  while 
in  an  incompressible  flow  (of  a  compressible  fluid),  the  effects  of  compress¬ 
ibility  may  be  assumed  to  be  equal  to  zero.  For  instance,  the  assumption  of 
incompressible  air  flow  is  acceptable  when  the  local  Mach  number  is  much 
less  than  one. 

Similarly,  the  flow  of  a  viscous  fluid  may  be  assumed  to  be  inviscid 
within  a  certain  region  when  effects  of  the  viscosity  are  negligible  in  that 
region.  Following  Prandtl,  the  flow  may  be  decomposed  into  a  boundary  layer 
(a  region  near  the  solid  surface,  i.e.,  the  surface  of  the  rotor  in  our  case)  where 
the  viscosity  is  important,  and  an  outer  region  where  the  effects  of  viscosity 
may  be  negligible:  we  will  concentrate  on  the  outer  region,  and  initially 
assume  that  the  flow  is  inviscid  in  that  region.  [However,  as  mentioned 
above,  the  issue  of  the  validity  of  the  inviscid-flow  assumption  in  the  region 
of  the  wake  may  be  too  restrictive  in  our  case:  this  point  is  addressed  at  the 
end  of  this  section,  where  we  will  use  the  ensemble-averaged  Navier-Stokes 
equations  to  discuss  how  the  formulation  may  be  modified  to  include  the 
presence  of  eddy-viscositv.] 

Under  the  above  assumptions  of  inviscid  incompressible  flow,  the  phe¬ 
nomenon  is  governed  by  the  continuity  equation  (conservation  of  mass)  for 
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incompressible  flows 


V  -v  =  0 


and  Euler  equations  (conservation  of  momentum) 


Dv 

Dt 


(3.1) 


(3.2) 


3.2  Helmholtz  Decomposition 

The  method  used  here  is  based  on  the  classical  Helmholtz  decomposition 
(see  Ref.  30)  of  the  flow  field  into  an  irrotational  part  and  a  solenoidal  part 

V  =  Vo  +  V  v  A  (3.3) 

In  Equation  3.3  <t>  is  the  scalar  potential  and  A  is  the  vector  potential,  which 
satisfies  the  condition 

V  •  A  =  0  (3.4) 

and  the  equation 

V2A  =  -«*  (3.5) 

where 

w  =  V  x  V  (3.6) 

is  the  vorticity. 

[In  order  to  prove  Eq.  3.1,  note  that  for  any  vector  » 


V  >  V  *  a  =  V(V  •«)  -  V2a  (3.7) 

Hence,  using  Eqs.  3.4  (which  may  always  be  satisfied  by  using  A  =  A0  *  V* 
with  V2x  =  -  V  •  A0  to  3.7.  one  obtains 

y  „  (v  -  V  ■  A)  =  w  -  T(V  •  A)  -  V2A  =  0  (3.8) 

which  indicates  that  there  exists  a  potential  function  <t>  such  that 

V  -  V  *  A  -  V<s  (3.9) 


in  agreement  with  Eq.  3.3.] 


Note  that  combining  Eqs. 
ible  flows, 


3.1  and  3.3,  one  obtains  that,  for  incompress- 
v2«,  =  0  (3.1°) 


3.3  Biot-Savart  Law 

If  w  =  0.  a  particular  solution  of  Eq.  3.5  is  (see  Eq.  2.51  applied  to  the 
infinite  space) 

13,11 

where 

r  -  v  -  x  (3-12) 


[Note  that 

v  A'fv'  v(i)Jl,(y'  *  -  HL~M^)dV[y) 

■  |315, 

where  n  is  the  outer  normal  to  a.  Noting  that 

V  -w  =  V- V  x  v  =  0  (3-14) 


we  can  see  that  Eq.  3.4  (used  in  Eq.  3.8)  is  not  satisfied  unless 

i#  •  n  =  0  on  <r 


(3.15) 


If  this  condition  is  not  satisfied,  then  as  shown  in  Ref.  30.  a  distribution  of 
vorticity  may  be  added  inside  a  to  generate  a  combined  solenoidal  vorticity 
field  in  the  whole  space  such  that  Eq.  3.4  is  automatically  satisfied,  the  last 
term  in  Eq.  3.13  does  not  exist  in  this  case.] 

It  should  be  noted  that  Eq.  3.11  is  fully  equivalent  to  the  classical  Biot- 
Savart  law.  In  order  to  show  this,  let  us  consider  a  very  thin  vortex  tube.  In 
this  case,  we  have  that  the  velocity  induced  by  the  vortex  tube  is  given  by 
(note  that  dV  =  dads,  where  a io  is  an  element  of  a  cross-section  of  the  vortex 
tube,  normal  to  w.  and  ds  is  an  arclength  parallel  to  •*.  so  that  uds  —  dy ) 


(3- 16) 
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Assuming  that  a  is  infinitesimally  small  (so  that  (y  -  x), >3  may  be  taken 
outside  of  the  surface  integral)  and  noting  that  the  intensity  of  the  vortex 
tube 

r  =  JJ  udo  (3  17) 

is  independent  of  s.  one  obtains  Biot-Savart  law: 

v(x)=r/c^Xiy  (318) 

Hence.  Eq.  3.11  may  be  considered  as  the  Biot-Savart  law  for  continuous 
distribution  of  vorticity. 


3.4  Method  of  Solution 

The  Helmholtz  decomposition  theorem  and  Biot-Savart  law  for  distribu¬ 
ted  vorticity  may  be  combined  to  yield 


V  =  V©  -r  W 


(3.20) 


where  w  is  the  velocity  induced  by  the  vorticity  and  is  given  by 

w(x)  =  v  *  IL  i^dV{y)  13211 

(with  «#  •  n  =  0  on  the  boundary  of  V.  see  Eq.  3.15).  whereas  ©  satisfies 
Laplace's  equation 

V2©  =  0  (3  22) 


The  condition  at  infinity  is  ©  =  0.  Therefore,  ©  satisfies  Eq.  2.59.  without 
the  wake  term  since  ©  is  continuous  in  the  whole  field: 


£(x)o(x) 


(3  23) 


Next,  a  boundary  condition  on  d<t>/dn  is  needed.  In  order  to  impose  this 
boundary  condition,  we  follow  the  classical  approach  of  Lighthill  (Ref.  37) 
and  introduce  a  "transpiration  velocity"  (i.e..  “equivalent  sources"  of  Ref. 
37):  this  velocity  is  such  that  an  inviscid  flow  around  a  permeable  surface 
(with  flow-through  defined  by  such  transpiration  velocity)  is  equal  to  the 
outer  flow  region  of  the  actual  viscous  flow  around  an  impermeable  surface. 


1 
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The  transpiration  velocity  for  three-dimensional  flows  was  introduced  by 
Lighthill  (Ref.  37)  and  is  given  by  Eq.  25  of  Ref.  37,  or.  using  our  notations. 


vtr  =  V,  • 


«(c)) dc 


(3.24) 


where  V,-  is.  the  surface  divergence.  6  is  the  boundary  layer  thickness,  and  < 
is  an  arclength  normal  to  the  boundary  surface. 

In  this  work,  we  will  assume  that  the  transpiration  velocity  is  known 
(using,  for  instance,  a  boundary-layer  code  or  a  thin-laver  Navier-Stokes 
code) .  The  boundary  condition  is  then  given  by  v  •  n  =  vb  ■  n  -  vtr .  or 


do 

—  =  -w  •  n  -  vb  ■  n  -  vtr 
On 


(3.25) 


[In  the  actual  numerical  calculations,  the  transpiration-velocity  contribution 
is  neglected:  in  this  case.  Eq.  3.25  is  equivalent  to  Eq.  2.6.] 

It  is  apparent  that  in  order  to  complete  the  formulation,  we  need  an 
equation  for  w  (as  well  as  suitable  initial  and  boundary  conditions).  The 
differential  equation  for  »  is  immediately  obtained  by  taking  the  curl  of  the 
Euler  equations.  By  using  the  identities 


(v  .  V)v 


<*  > :v 


and 


V  *  (a  <  b)  =  -(a  •  V)b^  (b-  V)a  -  a(V  •  b)  -  b(V  -  a) 
one  obtains  the  Helmholtz  vorticity  equation  (Ref.  30.  p.  152) 


(3.26) 


(3.271 


^  =  («-S>  (3.28) 

The  solution  of  Eq.  3.25  as  well  as  the  initial  and  boundary  conditions  are 
discussed  in  Sections  3.6  and  3.7. 


3.5  Material  Contravariant  Components 

The  objective  of  the  next  three  sections  is  to  introduce  the  concept  of 
material  (or  convected)  contravariant  components  of  vorticity  (that  is.  the 
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contravariant  components  of  the  vorticity  in  a  system  of  curvilinear  coordi¬ 
nates  that  are  convected  by  the  fluid  flow),  in  order  to  facilitate  the  integra¬ 
tion  of  the  vorticity  transport  equation.  In  particular,  it  will  be  shown  that, 
for  inviscid  incompressible  flows,  the  material  contravariant  components  of 
the  vorticity  are  constant  following  a  material  point.  This  result  is  used 
to  obtain  a  simple  but  powerful  computational  scheme  for  the  solution  of 
Helmholtz  vorticity  equation.  The  formulation  presented  here  is  based  on 
the  results  obtained  in  Refs.  35  and  36  for  compressible  flows. 

In  order  to  introduce  the  concept  of  material  contravariant  components, 
consider  a  system  of  material  (or  convected)  coordinates.  i.e..  a  system 
of  curvilinear  coordinates  that  are  convected  with  the  material  particles. 

As  in  the  classical  formulation,  a  given  material  particle  is  identified  by 
the  material  coordinates  £“  and,  at  any  time  its  location  is  determined  by 
the  Cartesian  coordinates  r,,  functions  of  the  material  coordinates  and  time 

A  =  z.Ua  0  (3  29) 

Equation  3.29  gives  the  Lagrangian  description  of  the  fluid  motion.  It 
may  by  noted  that  if  the  £“‘s  (a  =  1,2,3)  are  kept  constant,  then  z,  =  z,(t) 
represents  the  trajectory  of  a  fluid  point.  As  a  consequence,  coordinate 
lines  and  coordinate  surfaces  (for  instance,  the  surface  £'  =  const)  are  always 
composed  of  the  same  particles,  and  therefore  are  material  lines  and  material 
surfaces  respectively. 

It  should  be  emphasized  that  the  coordinates  are  closely  related  to 
the  classical  Lagrangian  coordinates  E“  (see  Ref.  30.  p.  128)  which  coincide 
with  the  Cartesian  coordinates  x,  at  t  =  0.  The  only  difference  between  the 
E“-coordinates  and  the  ^“-coordinates  is  that,  in  contrast  with  the  classi¬ 
cal  approach,  we  will  not  assume  that  the  ^“-coordinates  coincide  with  the 
Cartesian  coordinates  at  time  t  =  0.  This  not  only  emphasizes  the  curvilin¬ 
ear  nature  of  the  E“-coordinates  and  ^“-coordinates,  but.  more  importantly. 

*  Note  that  Greek  letters  are  used  for  curvilinear  components  (subscript  for  covariant  and 
superscript  for  contravariant  components).  Latin  letters  (subscript)  are  used  for  Cartesian 
components.  Einstein  summation  convention  on  repeated  indices  is  used  on  both  Greek  and 
Latin  indices. 
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yields  an  additional  flexibility  which  allows  for  a  convenient  choice  of  the 
^-coordinates  at  t  -  0  and  facilitates  the  derivation,  interpretation  and  ex¬ 
tension  of  some  useful  classical  results. 

Next,  some  elementary  concepts  on  curvilinear  coordinates  are  briefly 
reviewed  and  applied  to  the  specific  case  of  material  curvilinear  coordinates. 
The  velocity  of  a  material  point  is  given  by 


(where  after  a  vertical  bar  indicates  that  the  derivative  is  taken  with  (a 
=  constant).  Next  consider  a  quantity  /.  such  as  the  density  p.  which  is  a 
function  of  space  and  time.  The  material  (or  substantial)  derivative.  Df  Dt 
(see  Eq.  2.3),  is  the  time  derivative  for  an  observer  that  moves  with  the  fluid 
particle,  i.e.. 

Df  .  df  3f  df  dx. 


dx ,  dt 


or.  using  Eq.  3.30 


Df  df 

=  -t-  -  v  •  Vf 
Dt  dt 


Next,  for  any  given  time,  t ,  introduce  the  base  vectors 


which  are  tangent  to  the  coordinate  lines.  Recalling  that  the  £“-lines  are 
material  lines,  we  see  that  the  vectors  ga  are  always  tangent  to  the  same 
material  line.  Note  that 


D%a  _  d_  (  dx  \ 
Dt  dt  V3£a  / 


This  relationship  is  the  key  to  the  interpretation  of  the  vortex  stretching 
term  in  the  vorticity  evolution  equation. 

Next,  recall  that  any  vector,  in  particular  the  vorticity  w  may  be  ex¬ 
pressed  as 

«=w°gQ  (3  35) 

where  the  u.  a  's  are  known  as  contravariant  components  of  the  vorticity  vector 
«.  In  order  to  emphasize  that  the  £a's  are  material  curvilinear  coordinates, 
we  will  refer  to  the  /' s  as  material  contravariant  components  of  the  vorticity. 


wm  -  m  W  w  v »  „  W'J  .  *  IV  »  .  *  .  ■  V  ■  .  ■ — J  w  -  ^  .■■!-.  -  ■  r  •  J  - - p-  ■  i-  ^  v- 


3.6  Vorticity  Transport  in  Materia]  Contravariant  Components 


Combining  Eqs.  3.28  and  3.35,  we  obtain 


D_ 

Dt 


(w“g„) 


(3.36) 


The  second  term  in  Eq.  3.36  is  known  as  the  vortex-stretching  term: 
using  Eq.  3.34.  we  may  see  that  the  vortex-stretching  term  may  be  reinter¬ 
preted  as  representing  the  stretching  of  the  material  base  vectors.  Combining 
Eqs.  3.34  and  3.36.  one  obtains 


which  shows  that,  for  inviscid  incompressible  flows,  the  material  contravari¬ 
ant  components  of  the  vorticity  are  constant  following  a  material  point. 
Hence,  in  general,  the  vorticity  is  given  by 


-(f'5.*)  *.(€'•«)  (3  38) 

[If  £“  coincide  with  the  coordinate  H°,  (i.e..  coincides  with  the  Cartesian 
coordinates  x,  at  t  =  o)  then  Equation  3.38  coincides  with  Eq.  17.5  of  Ref. 
30.  p.  152:  a  result  obtained  by  Cauchy  in  1815.] 

The  above  results  are  used  here  as  the  basis  of  the  computational  method 
for  rotational  flows.  For  incompressible  inviscii  flows,  the  vorticity  is  given 
by  Eq.  3.38  at  all  times:  assume  that  we  know  the  vorticity  field  at  time  t  =  t(, 
and  the  functions  x(£a.z)  at  time  t0  and  at  time  ti ;  it  should  be  emphasized 
that  whereas  these  pieces  of  information  are  not  sufficient  to  evaluate  the 
velocity  field,  nontheless  they  are  sufficient  to  evaluate  the  vorticity  field 
at  time  t|.  Indeed,  from  x($“,t0)  we  may  evaluate  the  vectors  ga  at  t  =  <c, 
and  hence,  the  contravariant  components  *■“(£*.  to).  On  the  other  hand,  the 
knowledge  of  x<£°  t,)  allows  for  the  evaluation  of  ga  at  (  -  t,  and  hence,  of 
the  vorticity  field  at  time  t  =  tt,  via  Eq.  3.38. 

In  particular,  for  the  helicopter  rotor  wake  problem  of  interest  here,  once 
the  locations  of  materia!  gridpoints  on  the  wake  are  evaluated  at  a  given 
time  step,  the  direction  and  intensity  of  the  vorticity  field  can  be  evaluated 
immediately  using  Eq.  3.38.  Once  the  vorticity  is  known,  the  velocity  field 
is  evaluated  using  the  formulation  presented  in  Section  3.4. 


-  9 


It  may  be  worth  noting  that  some  classical  vortex  theorems  may  be 
obtained  (more  directly  than  in  the  classical  approach.  Ref.  30)  as  an  im¬ 
mediate  consequence  of  Eq.  3.38.  The  first  (and  most  relevant  here)  is  that, 
for  inviscid  incompressible  flows,  a  vortex  line  is  a  material  line.  In  order 
to  prove  this,  let  us  choose  the  coordinates  such  that  the  ^-lines  (and 
hence,  the  base  vectors  gi)  are  parallel  to  the  vorticitv  field  at  time  t  =  t0. 
This  implies  =  *3  =  0  at  i  =  t0  and  according  to  Eq.  3.37.  =  0  at  all 

times,  or 

-(^.t)  =  w1(^.to)  g.(^.«)  (3  39) 

Hence,  u  is  always  parallel  to  the  -lines,  which,  as  pointed  out  above, 
are  by  definition  material  lines:  therefore,  vortex  lines  are  material  lines  (in 
agreement  with  Ref.  30,  p.152). 

[Another  classical  result  that  is  easily  obtained  from  Eq.  3.38  is  related 
to  vortex  stretching.  Using  Eq.  3.39.  one  obtains 

---  =  k  3  40) 

fclcj 


where 

*  =  gi ,  gi  o  =  dx  !  dx  o  (3.41) 

is  the  stretching  factor  of  the  vortex  line  at  the  given  material  point.  Finally, 
it  may  be  worth  noting  that  if  we  introduce  a  small  material  cylinder  around 
an  element  dx  of  a  vortex  line,  having  a  cross-section  dA  and  a  constant  mass 
dm  -  pdA  dx  =  p<M0  dx ic.  one  obtains,  using  Eqs.  3.-41  and  3.42. 

udA-  uo'dAo  (3-42) 

in  agreement  with  Kelvin's  theorem,  (see  Eqs.  2.7  and  2.8). 

As  mentioned  above,  all  these  results  are  valid  for  compressible  isentropic 
flows  as  well  (see  Ref.  36).] 
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3.7  Thin  Wake  Approximation 

In  this  section,  we  will  show  how  the  formulation  can  be  simplified  in  the 
case  of  a  thin  wake.  i.e..  in  the  case  in  which  the  vorticity  in  the  outer  region 
may  be  assumed  to  be  approximately  equal  to  zero  everywhere,  except  for 
a  thin  (but  finite-thickness)  layer.  In  this  case,  the  velocity  induced  by  the 
vorticity  is  given  by  Eq.  3.21,  or 

w(x)  ’ v  1  ffL  Zr,iV[y]  ’  v  *  L„  ($.  *  15  451 


where  a  is  a  vortical  surface  (i.e.,  w  is  tangent  to  cr),  and  s-  is  an  arc  length 
in  the  direction  normal  to  a.  The  origin  of  f  is  “in  the  middle  of  the  layer.” 
and  is  defined  precisely  later  (see  Eq.  3.47).  The  thickness  6  is  such  that, 
in  any  case,  |w{  =  o  for  fi  >  6.  If  the  wake  is  thin,  then  Eq.  3.43  may  be 
approximated  with 


( 3.44 ) 


where 


7  = 


hide 


(3.45) 


and  a  is  the  midsurface  (f  =  o)  of  the  wake-layer. 

Next,  let  us  choose  the  material  coordinates  so  that  the  {Mines  are  par¬ 
allel  to  the  vorticity  field.  Then  u>2  =  ws  =  0  and  •#  is  given  by  Eq.  3.39. 
Also,  the  {Mines  are  chosen  so  that  the  surface  a  is  composed  of  a  network 
of  {'-lines  and  {Mines  (i.e.,  a  is  given  by  {s  =  o).  The  lines  {*  are  such  that 
at  time  t  -  t0.  they  are  normal  to  a.  Combining  Eq.  3.39  with  Eq.  3.45.  one 
obtains 


7  = 


(3.46) 


Next  let  us  choose  the  origin  of  o  so  that 


y<;d(  =  o 


(3.47) 


and  assume  that  g!  is  approximately  given  by  gj  =  at  -<■!>!  *  0({2)  (where  a,  is 
the  base  vector  gi  of  the  midsurface  a  (c  =  o).  and  b]  is  an  arbitrary  vector). 


I 
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Combining  with  Eqs.  3.46  and  3.47  and  neglecting  terms  of  order  S2.  one 
obtains 

7  =7**1  (3.48) 


where 


/6I2 

W1 

-6,2 


dc 


(3.49) 


Note  that  by  definition  of  the  wake  thickness.  =  0  at  the  extremes  of 
integration.  Therefore,  applying  Leibniz  rule 


It 


one  obtains,  using  Eq.  3.37, 


Dt1 

Dt 


dui1 

It 

dw* 

It' 


i  * ,  {  J  ,ti  —  contt 


i  1  .^2.^3=con«t 


(3.50) 


(3.51) 


[Note  that  the  ~  sign  is  due  to  the  fact  that  the  f-lines  are  always  normal  to 
a,  whereas  the  £S-Iines.  initally  normal  to  <r,  move  with  the  flow,  and  to  move 
away  from  the  normal  because  of  the  shear-motion  induced  by  the  vorticity. 
However,  this  approximation  used  in  Eq.  3.51  is  acceptable  if  du'  d(*  is 
small  (compared,  for  instance,  to  wl|m„/rf).] 

It  may  be  worth  noting  that,  at  least  in  the  case  of  a  flat  constant¬ 
thickness  layer  of  constant  vorticity,  Eq.  3.45  gives  the  exact  velocity  dis¬ 
tribution  outside  the  layer.  For  instance,  if  w  =  uj  (with  w  =  const )  for 
-6/2  <  f  <  b:2  and  x  =  0,  otherwise,  one  obtains  v  =  - 2  for  ;>  6/ 2.  (Actu¬ 

ally.  Eq.  3.45  gives  the  exact  result  even  if  u  =  *(c)j  with  w  =  0  for  <■  >  6/ -2  ). 
Therefore.  Eq.  3.45  is  expected  to  give  good  results  as  long  as  the  curvature 
is  not  too  high.  [This  is  not  the  case  near  at  high  tip  vortices.  However,  in 
that  case,  one  could  use  the  analysis  of  Section  3.3  (Eq.  3.16  and  3.18)  to 
reach  similar  conclusions  about  the  validity  of  Eq.  3.45.] 

Of  course,  the  velocity  obtained  from  Eq.  3.45  would  be  incorrect  inside 
the  wake  layer.  Therefore,  it  is  important  to  obtain  an  equation  that  gives 
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w,  in.  - 


the  time-evolution  of  the  wake  thickness.  Consider  the  second  moment  of 
vorticity 


7  =  /  r2"!?)^?  =  7**i 

J- 6,2 


|S  52) 


where 


,6/2 

71  =  /  f2"1  (?)<*? 

J -6/2 


(3.53) 


and  note  that  (with  the  same  qualifying  statement  used  for  Eq.  3.51) 


(3.54) 


In  order  to  obtain  an  approximate  equation  for  the  evolution  of  the  thickness, 
we  will  assume  that,  as  the  thickness  changes,  the  shape  of  the  function 
describing  -dU)  does  not  change,  i.e., 


(3.55) 


[This  expression  may  be  thought  of  as  the  first  term  of  a  series  expansion  for 
^•‘(f).]  For  instance,  one  could  use  /{«)  =  0  for  |u  >  l  and  /(u)  =  l  -  2u3  ~  u* 
[or  f(u)  =  (l  *  cos  ttu ) / 2.  or  even  flu)  =  lj  for  !«:  <  1.  In  this  case, 

-1  ,62  /  ,6/2 

-7  =  /  s-2/(2<r/d)<f<r  /  /  /(2f/<)if  =  «62  (3.56) 

~t‘  J-6  2  /  J -6/2 


where 


a 


u3  f[u)du 


(3.57) 


is  a  factor  independent  of  6  and  time;  if  /(«)  =  (1  -  2u2  -  u4),  then  a  =  1  7  [if 
f[u)  =  (l  -  cos  jtu)  2.  then  a  =  |  if  /(u)  =  1,  then  a  =  1/3].  Therefore  using 
Eqs.  3.51  and  3.54.  one  obtains 


that  is,  that  the  wake  thickness  remains  constant  in  time. 


(3.58) 
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3.8  Summary  of  Formulation  for  Inviscid  Rotational  Flows 

The  solution  of  Euler  equations  in  the  region  outside  the  boundary  layer 
has  been  reduced  to  the  solution  of  the  vorticity  transport  equation.  This 
may  be  studied  by  a  Lagrangian  analysis  of  the  points  of  the  vortical  region 
(wake),  i.e.,  by  solving  for  the  function  x  =  x(£“,t).  From  this  function,  one 
obtains  the  base  vectors  ga  =  dx/d(a  and  hence,  the  vorticity  w  =  u;Qg0,  where 
u'°  is  constant  in  time  and  equal  to  the  value  it  had  at  time  t  =  0  or  when  the 
fluid  particle  left  the  inner  region  (at  the  trailing  edge) .  From  the  distribution 
of  u.  one  may  evaluate  the  velocity,  w,  “induced’"  by  the  vorticity,  according 
to  Biot-Savart  law.  Then  a  potential  component,  V<s.  is  added  to  satisfy  the 
normal  boundary  condition.  This  process  yields  a  relationship 

^x(r.t)=vix(r,t),t;  (3.59) 

which  is  a  differential  equation  for  the  unknown  function  x(£°.t)-  This  equa¬ 
tion  in  general  is  solved  numerically,  using  for  instance  the  proced  ure  outlined 
in  Section  4. 

For  attached  flows,  one  may  use  a  thin-wake  approximation  in  which  only 
the  geometry  of  the  midsurface  of  the  wake  is  required.  This  implies  that, 
in  Eq.  3.59  only  two  curvilinear  coordinates  are  used  (i.e..  a  =  1.2  instead  of 
o  =  1,2,3).  Then  is  obtained  from  7  =  7“aa,  where  a„  are  the  surface  base 
vectors  (and  a  again  ranges  from  1  to  2).  The  velocity  outside  the  wake  layer 
is  obtained  using  Biot-Savart  law  for  a  layer  of  vorticity  (whereas  a  simple 
correction  is  required  inside  the  wake  layer,  which,  as  shown  above,  has  a 
thickness  that  remains  constant  in  time). 

3.9  Effect  of  Eddy-Viscosity  -  Navier-Stokes  Equations 

It  seems  natural  at  this  point  to  ask  ourselves  whether  the  assumption 
that  the  effects  of  the  viscosity  are  negligible  is  physically  acceptable.  More 
realistically,  one  should  speak  of  the  eddy-viscosity  since  the  wake  is  tur¬ 
bulent  and  typically  the  eddy-viscosity  coefficient  is  much  higher  than  the 
viscosity  coefficient.  In  order  to  discuss  this  issue,  it  is  convenient  to  see  how 
the  formulation  may  be  modified  by  using  the  ensemble-averaged  Navier- 
Stokes  equations  (instead  of  Euler  equations)  to  study  the  region  outside 
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the  boundary  layer  (i.e.,  without  the  restrictive  assumption  of  zero  eddy- 
viscosity). 

In  this  case.  Euler  equations.  Eq.  3.2.  are  replaced  by  the  ensemble- 
averaged  Navier-Stokes  equations* 

^  =  (3.60) 

where  uE  is  the  eddv-viscosity  coefficient.  Equations  3.3  and  3.24  are  still 
valid  whereas  Eq.  3.25  is  replaced  by  (assuming  for  simplicity  that  vE  is 
constant) 

^  =  (w  •  V)v  -  i/£V2«  (3.61) 

(A  less  restrictive  derivation  of  Eq.  3.61  (not  requiring  the  assumption 
that  j/£  is  constant)  may  be  obtained  directly  from  the  vorticity  evolution 
equation  for  incompressible  viscous  fluids  (i.e.,  Eq.  3.61  with  v  instead  of 
i/E).  It  should  be  emphasized  that  this  equation  indicates  that  »  is  not 
destroyed,  but  only  diffused  by  the  presence  of  viscosity.  Equation  3.61  is 
then  obtained  by  averaging  and  setting*  {*#'  •  Vv'  -  v'  •  V«')  -  i/V2«  =  i/EV2**: 
this  last  assumption  is  justified  by  the  fact  that  the  term  within  does 
not  destroy,  but  only  convects  (diffuses)  vorticity.] 

It  should  be  emphasized  that  the  presence  of  eddy-viscosity  does  not 
affect  the  method  except  for  the  equation  of  vorticity  evolution.  In  this  case, 
the  analysis  of  Section  3.6  must  be  generalized  to  viscous  flows.  In  particular. 
Eq.  3.36  becomes 

Dwa 

-^-g«  =  •'E  V2**  (3  62) 

or 

^  =  I'E  g°  ■  V2«  (3.63) 

where  ga  are  the  contravariant  base  vectors  such  that 

g°  •  g3  =  (3.64) 

where  <5j  is  the  Kronecker  delta  (an  extension  of  Eq.  3.63  to  compressible 
flows  is  given  in  Ref.  36). 

*  In  this  Section,  v  and  w  indicate  ensemble-averaged  values. 

Here  v'  and  indicate  perturbation  from  the  ensemble-averaged  values;  (  ;  indicates 


averaging. 


Equation  3.63  indicates  that  the  only  effect  of  the  viscosity  in  the  outer 
region  is  that  the  material  contravariant  components  of  the  vorticity  are  not 
simply  convected,  but  also  diffuse  in  the  process  of  being  convected. 

A  more  quantitative  result  may  be  obtained  for  attached  flows,  under  the 
assumption  of  thin  wake.  In  this  case,  we  may  assume  that  the  derivatives 
of  u  in  the  normal  direction  are  much  larger  than  those  in  the  tangential 
direction,  so  that  at  any  instant 

(3.65) 

Using  Eqs.  3.46,  3.63  and  3.65.  and  noting  that  (since  dm/d f  =  0  outside  the 
wake) 

f6/ 2  d2m  dm  6/2 

/  Tid<  =  S~  =  0  (366) 

J-6/ 2  -6/2 

one  obtains 

-W  =  0  (367) 

i.e..  that  Eq.  3.51  is  valid  in  presence  of  eddy-viscosity  as  well.  This  result 
should  have  been  expected  because:  (i)  the  effect  of  eddy-viscositv  is  not  to 
dissipate,  but  to  diffuse  the  vorticity,  and  (ii)  the  only  diffusion  included  in 
Eq.  3.65  is  that  in  the  direction  which  does  not  affect  t  (the  integral  of  m 
in  the  c  direction,  Eq.  3.45). 

For  the  same  reason,  one  should  expect  that  the  time  evolution  of  the 
thickness  of  the  wake  is  affected  by  the  presence  of  eddy-viscosity.  In  order 
to  show  this,  multiply  Eq.  3.63  by  c2  and  integrate  to  obtain  (choosing  the 
{Mines  parallel  to  w.  using  Eq.  3.53  and  neglecting  higher  order  terms) 


D71 ,  f 

~dt  1  =  UE  J_ 


2d2u 

d? 


Noting  the 
3.48) 


d$  -  0  outside  the  wake,  one  obtains  (see  Eqs.  3.45  and 


Z"4'2  d2m  f 4/2  dm  f4’2 

/  =  -  /  2<^*  =  2  /  =  27  =  2V«, 

J-6  2  d<2  J-6! 2  3?  J- 6/2 


and  hence,  combining  with  Eq.  3.68. 


=  2  I 'El1 
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Finally,  using  Eqs.  3.56  and  3.67,  and  recalling  that  Da/dt  =  0  (see  Eq.  3.57), 
one  obtains 

=  a-f1  ^-  =  2  i/£  Tf1  (3.71) 

or.  assuming  for  simplicity  that  uE  is  constant, 

SU'.t)  =  (*g  +  2  uE  </q)5  (3.72) 


with 

*>  =  ftf'.O)  (3.73) 

where  t  =  o  is  either  the  initial  time  or  the  time  when  the  wake  point  left  the 
inner  region  at  the  trailing  edge. 

The  approximate  formulation  presented  above  indicates  that  (within  the 
limits  of  the  formulation)  the  only  effect  of  eddy-viscosity  is  a  growth  of  the 
wake  thickness.  If  ReE  is  the  eddy-Reynolds  number  based  upon  the  eddy 
viscosity  coefficient,  chord  and  tip  speed  0 R 

cQR 


RtE  = 


vE 


one  obtains 


S0  /  a  \  S0 1  e  ReE 


(3.74) 


(3.75) 


If  So  is  of  the  order  of  the  thickness  of  the  blade,  or  c/S0  ~  10.  using  a  =  l  7 
(see  Eq.  3.57),  one  obtains  that  after  one  revolution  (i.e..  for  /  =  T  =  2v  n) 
the  thickness  increases  bv  the  factor 


£-^-  2800,5  ^ 

For  an  eddy  Reynold’s  number  of  5,000  and  R/c  -  16,  we  obtain  that 

6, 


So 


=  v'l  +  28.15  =  5.4 


(3.76) 


(3.77) 


i.e..  the  thickness  becomes  more  than  five  times  larger  during  the  first  revo¬ 
lution. 

In  conclusion,  the  effect  of  the  viscosity  may  be  taken  approximately  into 
account  by  introducing  a  thickness  growth  that,  while  not  extremely  small, 
is  sufficiently  small  to  be  evaluated  using  the  simplified  analysis  presented 
in  this  section.  It  should  be  emphasized  that  the  main  effect  of  the  wake 
thickness  is  “local.”  The  effect  upon  distant  points  may  still  be  evaluated 
using  the  zero-thickness  approximation  (see  Eqs.  3.43  and  3.44). 
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SECTION  4 


COMPARISON  OF  THE  TWO  FORMULATIONS: 
THE  LIMITING  CASE 


The  formulations  for  potential  flows  and  for  rotational  flows  have  been 
presented  in  Sections  2  and  3,  respectively.  The  first  formulation  is  limited 
to  zero-thickness  wake,  whereas  the  second  one  allows  for  finite-thickness 
wakes  and  may  be  used  for  both  inviscid  (Sections  3.1  to  3.8)  and  viscous 
flows  (Section  3.9).  It  is  apparent  that  the  limiting  case  of  the  rotational- 
flow  formulation  as  the  wake  thickness  goes  to  zero  (e.g.,  for  inviscid  fluids) 
should  be  fully  equivalent  to  the  potential-flow  formulation.  In  this  section, 
the  equivalence  between  the  two  formulations  (as  the  wake  thickness  goes  to 
zero)  is  analized  in  detail.  This  analysis  is  not  only  useful  for  the  theoret¬ 
ical  clarification  of  the  relationship  between  potential  flows  and  rotational 
(inviscid  and  viscous)  flows,  but  it  also  helps  in  indicating  how  the  numer¬ 
ical  scheme  used  for  the  potential-flow  formulation  may  be  modified  to  be 
applicable  to  the  rotational-flow  formulaltion. 

4.1  Rotational  Flows:  Limiting  Case 

As  the  wake  thickness  goes  to  zero  (i.e.,  the  viscosity,  and  hence,  the 
boundary-layer  thickness  go  to  zero,  see  Section  2)  the  formulation  for  rota¬ 
tional  flows  simplifies  considerably.  The  velocity  is  given  by  (see  Eqs.  3.20 
and  3.44) 


v  =  V<t>  +  w 


(4.1) 


with 


w(x)  =  V  x 


L. 


(4.2) 
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The  scalar-potential  <?  satisfies  the  equation  (Eq.  3.23) 


E(x)<t>(x) 


do  1 

dn  Akt  dn 


j 

J 


(4.3) 


with  boundary  conditions  (see  Eq.  3.25) 

do  .  . 

—  =Vi,n-wn  on  erb  (4.4) 

dn 

The  surface  distribution  of  vortiticy  7,  satisfies  the  equation  (see  Eq. 
3.52) 

7°  =  constant  in  time  (4.5) 

following  a  wake  point. 

Because  of  the  well  known  equivalence  between  doublet  layers  and  vortex 
layers  (Ref.  31),  Eq.  4.2  may  be  rewritten  as 


=  v/..a*4U;) 


(4.6) 


with  A <bw  such  that 

7  =  n  >  Vt(A®u,) 

(where  V£  is  the  tangential  gradient)  and 


(4.7) 


(4.8) 


Therefore,  the  rotational-flow  formulation  (in  the  limiting  case  of  zero-wake 
thickness)  may  be  restated  as  follows:  the  velocity  is  given  by 


v  =  V®  -r  V0W 


(4.9) 


where  ®  satisfies  Eq.  4.3,  whereas  ou  is  given  by  Eq.  4.8. 

Using  Bernoulli’s  theorem  and  noting  that  A®  =  0.  one  obtains  (see  Sec¬ 
tion  2.4) 

^(A®«,)  =  0  (4.10) 

or 

A®w  =  constant  in  time  (4.11) 


following  a  wake  point. 

[Note  that  Eqs.  4.5  and  4.11  are  equivalent.  In  order  to  show  this,  recall 
that,  in  general. 

V/  =  g»il-g2^.-g^  (4.12) 

1  8  dei  8  dp  8  '  ' 

where  g“  are  the  contravariant  base  vectors  which  are  such  that 

g"-g a  =  K  (a,  3  =1.2.3)  (4.13) 


This  implies 


g  ‘  g  =  g1  /  (gl  *  gs  •  gs) 


g  ”  g1  =  g  j  (gi  *  g2  -  gs)  (4.14) 

In  particular,  (using  «i  =  gi,  *2  =  g2,  =  a3  =  n  =  at  x  a2/;«i  x  a2  )  one  obtains 


7  =  n  x  VtfA^u,) 


=  V"Xa dV+n**‘ df*)**" 

(  d(Ad>u.)  d(A<M  \  / 


*2  )  /  *1  x  *2 


=  71*!  -  72a2 


where 


j  _  -1  d(Aou 

7  i«i  x  «2i  d£2 


I2  =  - - -  (4.16) 

Since  Aow  is  constant  in  time  following  a  wake  point  (i.e..  keeping  £*  and 
constant)  one  can  see  that  the  same  is  true  for  the  material  contravariant 
components  if  and  only  if 

>i  x  a2  =  constant  in  time  (4 . 17) 


following  a  wake  particle.  In  order  to  show  that  Eq.  4.17  is  correct,  consider 
a  small  material  wake-volume  element  dV  =  do6.  where  <*>  is  the  thickness  of 
the  wake  and  do  is  the  area  of  a  midsurface  element.  Since  the  densitv  is 


constant  and  dV  is  a  material  volume,  then  dV  is  constant  in  time.  But,  as 
shown  in  Section  3.  S  is  constant  in  time  for  inviscid  flows.  Therefore, 

da  =  ai  x  a2i  d^d^2  (4.18) 

is  also  constant  in  time,  and  hence,  Eq.  4.17  is  true  (since  d£l  and  d£2  are 
constant  for  a  material  element). 

It  may  be  worth  noting  that  from  a  physical  point  of  view.  Eq.  4.17  may 
be  stated  as  saying  that  the  area  of  a  material  element  of  the  wake  surface 
is  constant  in  time. 

In  conclusion,  the  material  contravariant  components  and  Ao„  are 
related  through  Eq.  4.16  and  therefore,  because  of  Eq.  4.17.  Eqs.  4.5  and 
4.11  are  equivalent.) 


4.2  Relationship  with  Potential-Flow  Formulation 


Presented  in  this  section  is  the  relationship  between:  (i)  the  limiting  case 
of  the  rotational-flow  formulation  (Eqs.  4.3  to  4.11),  and  (ii)  the  potential- 
flow  formulation  of  Section  2. 

In  order  to  facilitate  the  discussion,  consider  Eq.  2.59  and  note  that, 
following  Lamb  (Ref.  32),  this  equation  may  be  generalized  by  allowing  a 
'‘fictitious”  flow  inside  ab:  this  interiaor  flow  has  a  velocity  v;  =  Vy?;,  where 
ri  has  an  integral  representation  of  the  type 


Ei'Pi 


-t 


dni  4irr  dni 


-1 

4rr 


da 


(4.19) 


This  equation  is  obtained  by  following  the  same  procedure  used  to  obtain 
Eq.  2.52,  but  using  the  volume  inside  ab  instead  of  that  outside.  As  a 
consequence, 

E,  =  1  -  E  (4.20) 

and 

n  /  =  -n  (4.21) 

Combining  Eqs.  2.59  and  4.19  and  using  Eqs.  4.20  and  4.21.  one  obtains 

l 

4  srr 
d 


E<p  +  (1  —  E )'£>i  — 


e  / 1 


-(-) 

dn  \  4 nr  / 


da 


i?) 


da 


(4.22) 
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which  is  the  generalized  integral  equation  for  potential  flows.  Note  that  pj 
satisfies  Eq.  4.19  and  is  otherwise  completely  arbitrary  (since  Ej  =  0  outside 

<Th)- 

On  the  other  hand,  Eqs.  4.3  and  4.9  may  be  combined  to  yield 


Eo  —  <*>u. 


d0  1  _  d  (  Ml  _  ff  \  L 

dn  47r r  dn  \  4irr )  JJaa  u  dn 


(4.23) 


Note  that  (see  Eqs.  2.33  and  4.8)  the  velocity  is  given  by 


v  =  V<£  =  + 


(4.24) 


with  (see  Eqs.  2.36,  4.4  and  4.6) 

dip  d<t>  dip  w 
dn  dn  dn  k 

Therefore,  one  expects  that,  at  least  outside  a 


(4.25) 


p>  =  <p  <p  u, 


(4.26) 


Substituting  this  equation  into  Eq.  4.23,  one  obtains  (note  that  A©  =  o  on 
the  wake) 


Ep  -  (1  -  £)<?u 


p  -  a>w)- - [p  -  ®u)^- 

4  7r  r  dn 


da 


(4-27) 


By  comparing  Eqs.  4.22  and  4.27.  we  may  conclude  that  the  limiting  case  of 
the  rotational-flow  formulation  coincides  with  the  potential-flow  formulation 
as  extended  by  adding  an  internal  flow 


pl  =  0w 


(4.28) 


4.3  Comments 


The  results  of  Section  4.3  indicate  that  Eq.  2.59  (potential-flow  for¬ 
mulation)  and  Eqs.  4.3.  4.4  and  4.6  (limiting  case  of  the  rotational-flow 
formulations)  are  fully  equivalent.  However,  the  second  formulation  has  a 


distinct  advantage  over  the  first  one.  In  the  potential-flow  formulation,  the 
wake  is  treated  as  a  doublet  layer,  and  hence,  is  restricted  to  zero  thickness. 
In  the  limiting  case  of  the  rotational  formulation,  on  the  other  hand,  there 
is  no  way  to  ascertain  w'hether  the  wake  is  treated  as  a  doublet  layer  or  as  a 
vortex  layer.  This  implies  that  the  wake  thickness  does  not  have  to  be  equal 
to  zero,  and  therefore,  the  formulation  can  be  used  to  introduce  (in  a  gradual 
fashion!)  the  effects  of  viscosity  and  eddy-viscosity. 

It  may  be  worth  pointing  out  that  "artificial  viscosity"  may  be  introduced 
in  the  potential-flow  formulation  for  the  evaluation  of  the  velocity  of  the 
wake  points  (via  Eq.  2.69).  As  we  will  see  in  Section  5.  this  yields  results 
that  are  practically  indistinguishable  from  those  obtained  from  the  use  of 
actual  eddy-viscosity  in  the  rotational  flow:  however,  there  is  a  big  difference 
between  the  two  cases.  In  the  potential-flow  formulation,  “artificial  viscosity” 
is  just  what  it  says:  it  is  an  artificial  way  to  enhance  the  numerical  stability 
by  increasing  the  vortex  size  in  the  evaluation  of  the  velocity  of  the  wake 
points:  this  introduces  an  inconsistency  in  the  way  the  wake  is  treated  in  the 
integral  equation,  Eq.  2.59  (zero-thickness  doublet  layer)  and  in  the  velocity 
calculation.  Eq.  2.69  (vortex  layer  with  a  finite-thickness  stemming  from  the 
artificial-viscosity  approximation).  The  rotational-flow  formulation,  on  the 
other  hand,  is  fully  consistent:  the  same  (finite-thickness  or  zero-thickness) 
wake  description  is  used  for  both  the  integral  equation  (Eq.  3.22)  and  in 
the  velocity  expression  (Eq.  3.20  or  Eq.  4.2  for  finite-thickness  and  zero¬ 
thickness.  respectively).  Therefore,  the  effect  of  viscosity  is  not  necessarlv 
artificial,  but  may  be  the  one  actually  provided  by  the  turbulence  model  that 
is  chosen  to  represent  the  phenomenon. 

4.4  Discretization  of  Potential-Flow  Formulation 

In  order  to  solve  the  problem,  the  integral  formulations  presented  above 
must  be  discretized  in  space  and  time.  In  this  section,  we  consider  the  dis¬ 
cretization  of  the  potential-flow  formulation.  The  results  of  this  section  will 
then  be  modified  to  obtain  the  discretization  of  the  rotational-flow  formu¬ 
lation:  in  order  to  be  able  to  motivate  this  modification,  it  is  convenient 
to  introduce  a  finite-element  discretization  process  of  arbitrary  order  (even 


though  a  zeroth-order  formulation  is  used  in  the  actual  computations). 

Using  general  finite-element  representation,  with  M  nodes  on  the  sur¬ 
face  of  the  body  and  ,V  nodes  on  the  surface  of  the  wake,  it  is  possible  to 
approximate  p. 


(4  29) 


and  Ar-  as 


,w 

^(*'0  =  H  «(t)  -W,(x) 

t=l 


(4.30) 


M 

V»(x,t)  =  ^2  V*(«)  W.(x) 
1  =  1 


(4.31) 


and 


JV 

±p(x,t)  =  ^  Av?;(t)  Nj(x)  (4  32) 

j  =  i 


where  p,(t)  =  ^(x,,t).  *-,(t)  =  v>(x».<)  and  A*>,(<)  =  whereas  ,V/,(x)  and 

.Vj(x)  are  prescribed  finite-element  shape  functions/  which  have  the  property 


.Vf,  (x)  =  6,k  =  1  for  i  =  * 

=  0  for  t  *  h 


(where  6lk  is  the  Kronecker  delta).  For  instance,  in  the  zeroth-order  formu¬ 
lation  used  in  Refs.  1-3.  M,  and  .V,  are  equal  to  one  within  an  element  of  the 
surface  (of  the  body  and  the  wake)  and  equal  to  zero  outside.  For  a  first- 
order  formulation,  within  a  given  quadrilateral  element,  p  may  be  described 
as  (see.  for  instance.  Ref.  26) 

P=P—\[l  +  fl)d  *  fV 

P---M  +  fM(i  -  f2)- 
4 

f-.-M  -  f')(i  +  f2)- 

4 

*>-  Id-^Hl-f2)  (4  33) 

4 

*  Note  that  in  order  to  avoid  proliferation  of  symbols,  the  same  shape  functions  are  used 
for  p  and  t.  although  this  is  not  essential  to  the  method. 
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where  *c>+_,  <p~~ ,  £>__  and  <p--  are  the  values  of  -p  at  the  four  corners  of  the 
element.  (The  geometry  of  the  element  itself  may  be  described  in  a  similar 
fashion.)  Hence,  A/-_,  and  A/__,  the  local  shape  function  (which, 

assembled,  yield  the  global  shape  function  are  given  by 

M.-htf1.*2)  =  Jll  +  C1)!!  -£2)  (4-34) 

4 


and  similar  expressions  for  the  others. 

Combining  Eqs.  4.24  to  4.32  with  Eq.  2.59,  using  the  colocation  method 
and  setting  the  control  points  at  the  nodes  of  the  elements,  x*,  one  obtains 
(see  e.g.,  Ref.  26) 

N  N  M 

Ek<Pk  =  Bkii-' ,  +  Ckipi  -r  ^  Fkj  A(£j 

«  =  1  «=1  ]- 1 

where  Ek  =  E(xk)  and 

®“  "*(»)  <TF^nto,,i 

and 

-  J[/'w  K  (■ 

At  each  time  step,  v>,  are  known  from  the  boundary  conditions.  A*?;  are 
known  from  the  preceding  time  step,  and  Eq.  4.35  may  be  solved  for  y?*. 

A  similar  discretization  is  used  for  Eq.  2.69  to  obtain  the  velocity  of  the 
wake  points  in  terms  of  i\,  <px  and  A*?,.  Once  the  velocity  v;(t)  is  known,  the 
new  location  of  a  wake  point  is  obtained  by 

x(t  +  At)  =  x,(t)  *  Vj (r ) Ar  (4.39) 

Note  that  according  to  Eq.  2.38,  A is  constant  in  time  (as  long  as  j 
indicates  the  same  material  wake  point).  The  trailing  edge  condition  is  used 
to  obtain  the  values  of  Ay?,-  on  the  first  row  of  elements  from  the  values  of 
<px  at  the  trailing-edge  points  on  the  upper  and  lower  surfaces.  The  process 
may  then  be  repeated  for  the  next  time  step. 


(4.35) 

(4.36) 

(4.37) 

(4.38) 


For  the  zeroth-order  formulation,  the  numerical  algorithm  is  essentially 
identical  to  that  of  Refs.  1  to  3  (see  these  references  for  details).  The  main 
difference  between  the  present  formulation  and  that  of  Refs.  1  to  3  is  that 
here  the  first  wake  element  is  rigidly  connected  with  the  rotor  blade.  This 
has  two  advantages',  (i)  it  eliminates  the  unrealistic  trend  obtained  for  the 
first  element  in  Refs.  1  to  3,  and  (ii)  it  does  nor  require  a  matrix  inversion 
at  each  time  step  (as  needed  in  Refs.  1  to  3)  because  all  the  coefficients  of 
the  matrix  that  multiply  are  time  independent. 

4.5  Discretization  of  Rotational-Flow  Formulation 

The  discretization  for  the  rotational-flow  formulation  differs  from  that 
for  potential-flow  formulation  only  in  the  discretization  of  the  vorticity  (the 
integral  equation  may  be  discretized  in  the  same  way  for  both  problems, 
yielding  for  Eq.  4.3  an  approximation  equal  to  Eq.  4.35,  with  Fk,  -  0).  One 
of  the  difficulties  in  discretizing  the  vorticity  equation  is  to  impose  that  the 
vorticity  field  is  solenoidal.  This  problem  does  not  arise  in  the  potential- 
flow  formulation  because  an  arbitrary  doublet  distribution  corresponds  to  a 
vorticity  distribution  that  is  always  solenoidal  (Ref.  31).  The  relationship 
between  the  two  formulations  may  be  exploited  to  extend  the  discretization 
for  the  potential-flow  formulation  to  that  for  rotational  flows. 

In  order  to  accomplish  this,  it  is  convenient  to  start  from  a  zero-thickness 
formulation.  In  this  case  (see  Eqs.  4.7  and  4.32)  the  vorticity  may  be  ex¬ 
pressed  as 

7  =  7>  (4.40) 

» 

where 

7,=nxVN,(x)  (4.41) 

is  the  surface  distribution  of  vorticity  corresponding  to  a  doublet  distribution 
Ni(x).  Note  that  if  /V,  is  continuous. 7,  will  not  include  concentrated  vorticies. 
This  is  true,  for  instance,  for  the  first-order  formulation  which  yields  (see  Eqs. 
4.15  and  4.16) 

7-  -  =  -  |(1  -  f')»l  /  «1  ‘  *2 
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This  is  an  elementary  "vortex  sheet*"  with  vortex  lines  parallel  to  the  lines 
=  constant.  For  the  zeroth-order  formulation,  on  the  other  hand.  Eq. 
4.41  yields  a  quadrilateral  vortex  filament. 

In  order  to  extend  this  formulation  to  distributed  vorticity  (limiting  our¬ 
selves.  for  the  sake  of  simplicity,  to  the  case  of  a  thin  wake)  we  may  use  the 
concepts  introduced  in  Section  3.7,  Eq.  3.55,  and  set 

»-+(£l,£2,f)  =  jF( 2f/6)  nxVtNU\?)  (4.43) 

(with  f*  F(u)du  =  l)  which  also  automatically  satisfies  the  solenoidality  con¬ 
dition  (the  function  F(u)  may  be  obtained  from  the  expression  for  f(u)  given 
in  Section  3.7).  Equation  4.43  represents  a  thin  “blob”  of  vorticity  which 
automatically  satisfies  the  solenoidality  condition  and  automatically  adjusts 
the  intensity  to  take  into  account  the  vortex  stretching.  For  the  zeroth-order 
formulation,  we  are  back  to  the  expressions  for  finite-core  vortex  filaments 
(see,  for  instance,  Ref.  l). 


SECTION  5 

NUMERICAL  RESULTS 

In  order  to  assess  their  validity  and  efficiency,  the  two  formulations  pre¬ 
sented  in  Sections  2  and  3  have  been  coded  using  the  discretization  outlined 
in  Section  4.  As  mentioned  before,  only  the  zeroth-order  formulation  (i.e.. 
concentrated  vortices  or  equivalent  constant-strength  doublet  panels)  has 
been  implemented.  Numerical  results  have  been  obtained.  Detailed  presen¬ 
tation  of  the  results  is  given  in  Ref.  38,  which  is  meant  to  serve  as  a  collection 
of  results  and  will  be  continuously  updated  as  new  results  become  available. 
Only  those  results  that  are  needed  for  the  assessment  of  the  formulations  are 
included  here. 

All  the  results  presented  here  are  for  the  test  case  considered  by  Morino. 
Kaprielian  and  Sipcic  (Refs.  1  to  3)  and  by  Rao  and  Schatzle  (Ref.  7)  for 
which  the  experimental  results  of  Bartsch  (Ref.  16)  are  available.  The  main 
reason  for  this  choice  is  the  fact  that  we  are  familiar  with  this  problem  and 
prefered  to  work  on  familiar  grounds  while  developing  a  new  formulation 
and  a  new  computer  program.  It  is  understood  that  additional  results  on 
different  test  cases  are  needed.  They  are  currently  being  obtained  and  will 
be  included  in  Ref.  38.  However,  even  if  only  one  test  case  is  presented  here, 
the  results  presented  are  quite  extensive  because  of  the  many  issues  which 
had  to  be  addressed. 

5.1  Preliminary  Remarks 

The  test  case  presented  here  is  a  rotor  with  tip  radius  R  =  17.5  ft.  root 
cut-out  radius  r0  =  2.33  ft.  chord  e  =  1.083  ft.  collective  pitch  angle  6r  =  10.61 c 
and  twist  angle  =  -5C.  The  angular  velocity  is  n  =  355  rpm. 

Some  preliminary  studies  were  made  (using  a  single-bladed  rotor  with  a 
prescribed  wake)  in  order  to  assess  the  rate  of  convergence.  These  results  are 


presented  in  Ref.  38.  The  main  conclusion  is  that  the  two  formulations  yield 
very  similar  results  although  the  rate  of  convergence  of  the  rotational-flow 
formulation  is  not  as  fast  as  that  of  the  potential-flow  formulation  in  terms 
of  the  number  of  elements  on  the  blade.  In  any  event,  sufficiently  accurate 
results  are  obtained  for  both  formulations  using  7  elements  in  the  radial 
direction,  3  elements  in  the  chordwise  direction  and  a  five-spiral  wake  with 
12  elements  per  spiral.  All  the  results  presented  in  this  Section  use  the  above 
values.  Because  of  its  broader  applicability  (i.e.,  from  potential  to  viscous 
flows)  the  emphasis  here  is  on  the  rotational-flow  formulation;  the  potential- 
flow  formulation  is  only  used  to  assess  the  accuracy  and  the  efficiency  of  the 
rotational-flow  formulation. 

The  results  of  Ref.  1  show  that,  because  of  the  wake  truncation  the 
last  few  wake  spirals  tend  to  move  outward  and  upward  due  to  the  lack  of 
the  suction  effect  caused  by  the  portion  of  the  wake  that  was  ignored  in 
the  computational  model.  As  already  pointed  out,  for  instance  by  Summa 
(Ref.  20).  this  problem  may  be  corrected  by  introducing  an  intermediate 
and  a  far-wake  model.  Extensive  preliminary  studies  (presented  in  Ref.  38) 
were  made  in  order  to  assess  several  intermediate  and  far-wake  models.  The 
restrictions  that  we  imposed  on  the  selection  of  the  model  were:  (i)  the  results 
should  be  reasonably  insensitive  to  the  far-wake  model  employed  and  (ii)  the 
intermediate  wake  model  should  be  eventually  developed  directly  from  the 
computational  results  themselves. 

In  trying  to  develop  an  intermediate  wake  model,  we  were  not  able  to 
infer  the  shape  of  the  wake  from  the  computer  data  as  they  were  generated. 
However,  we  noticed  that  the  velocity  of  points  of  the  wake  were  quite  pre¬ 
dictably  regular  and  could  be  inferred  from  the  velocity  of  the  wake  points 
of  the  preceding  spirals.  Therefore,  we  chose  an  intermediate  wake  model 
in  which  the  velocity  (not  the  positions)  are  assumed.  Note  that  the  only 
function  of  the  intermediate  wake  spirals  is  to  enable  the  last  few  spirals  of 
the  wake  to  move  with  the  physically  expected  behaviour:  thus,  the  results 
obtained  for  the  sectional  lift  distribution  are  quite  insensitive  to  the  velocity 
assumed  in  the  intermediate  wake  spirals  because  the  exact  geometry  of  the 
last  few  spirals  is  not  as  important  as  the  preventing  of  their  divergence: 
aiso.  only  the  roll-up  in  the  first  one  or  two  spirals  is  crucial  in  predicting 


5-2 


the  blade  loading.  For  the  results  presented  here  (except  for  the  rotor  in 
ground-effect)  the  velocity  assumed  for  the  intermediate  spirals  is  the  one 
that  would  generate  a  Landgrebe  wake,  even  for  the  single-bladed  rotor  for 
which  the  Landgrebe  model  is  not  necessarily  valid.  (For  the  rotor  in  ground- 
effect.  the  velocity  for  the  intermediate  wake  is  described  in  Section  5.4.)  It 
is  acknowledged  that  an  approach  less  dependent  on  experimental  data  is 
desirable,  and  eventually  a  self-correcting  model  (i.e.,  one  with  parameters 
that  are  updated  during  the  running  of  the  code)  should  be  developed. 

The  far-wake  model  for  the  rotor  out  of  ground-effect  is  needed  to  com¬ 
pensate  for  the  wake  truncation.  It  is  known  (see  for  instance  Summa.  Ref. 
20)  that  a  uniform  vorticity  distribution  over  a  semi-infinite  cylinder  induces 
(outside  of  the  cylinder)  exactly  the  same  velocity  as  a  sink  disk  located  at 
the  end  of  the  cylinder,  provided  the  flux  into  the  sink  is  equal  to  the  flux 
through  a  cross-section  of  the  cylinder.  Since  the  flux  through  the  far  wake 
is  equal  to  the  flux  through  the  rotor,  the  total  flux  into  the  sink  must  be 
equal  to  trft2 «0.  where  u0  is  the  velocity  at  the  rotor  disk.  An  estimate  for  u0 
can  be  obtained  from  the  actuator  disk  theory  (for  a  rotor  in  hover)  accord¬ 
ing  to  which  u0  =  QR^Ct/2.  This  point  is  discussed  further  in  Section  5.2 
where  numerical  results  obtained  using  different  estimates  for  the  flux  are 
presented. 

Next,  consider  the  issue  of  the  initial  geometry  of  the  wake.  In  order  to 
assess  the  robustness  of  the  code,  the  initial  wake  was  assumed  to  have  little 
resemblance  to  the  expected  final  result,  i.e.,  a  modified  classical  wake  with 
a  pitch  As  =  2*RV'Ct,2  and  C>  prescribed  (see  the  following  subsections  for 
numerical  values).  The  modifications  are  outlined  in  the  following.  After 
some  experimentation,  it  was  observed  at  the  first  time-step  that  the  tip 
vortex  would  move  upward  (as  it  should),  but  the  vortex  filament  nearest  to 
the  tip  vortex  would  move  downward  and  take  a  long  time  to  recover  from 
this  position.  Therefore,  in  all  the  results  presented  here,  the  initial  classical 
wake  was  modified  by  moving  the  tip  vortex  upward  by  one-fourth  of  the 
pitch  (but  not  higher  than  the  blade  tip).  Moreover,  a  radial  contraction  of 
the  wake  was  employed  (using  Landgrebe  s  expressions  even  for  the  case  of 
a  single-bladed  rotor  for  which  they  are  not  strictly  applicable.) 
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In  addition  to  the  initial  geometry  of  the  wake,  there  is  the  issue  of 
the  initial  vorticity  distribution  on  the  wake.  In  order  to  avoid  the  well- 
known  problems  of  the  start-up  vortex  (see.  e.g..  Summa,  Ref.  6)  in  all  our 
studies,  computations  for  the  first  .V  time-steps  were  performed  with  the 
wake  frozen  in  the  initial  position.  (Here  .V  is  obtained  from  the  number  of 
elements  along  a  vortex  filament,  by  adding  two,  e.g..  ,V  =  62  for  five  spirals 
and  twelve  elements  per  spiral.)  This  enables  ‘‘loading”  of  the  wake  w'ith 
a  vorticity  distribution  very  close  to  that  of  the  frozen-wake  steady-state 
solution  and  ensures  that  any  start-up  vorticity  spike  is  washed-out  from  the 
wake. 

Another  issue  that  we  addressed  in  the  preliminary  studies,  not  included 
here  (see  Ref.  38),  is  that  of  the  vortex  at  the  root  cut-out.  The  vortex 
sheet  tends  to  roll-up  at  the  root  just  as  it  does  at  the  tip.  It  is  necessary  to 
employ  small  elements  at  the  root  regie  a  in  order  to  calculate  this  roll-up. 
However,  it  is  well-known  (see  e.g.,  Miller.  Ref.  21)  that  this  region  is  not 
crucial  to  the  overall  phenomenon.  Therefore,  we  used  large  elements  near 
the  root.  As  a  consequence,  the  vortex  sheet  does  not  roll-up.  but  tends 
to  move  upward  and  generate  disturbances  that  eventually  propagate  to  the 
entire  wake.  In  order  to  avoid  this,  the  first  vortex  of  the  wake  sheet  has 
a  prescribed  velocity  (i.e..  we  use  the  same  methodology  employed  for  the 
intermediate  wake.) 

Before  discussion  of  the  results,  a  few  words  about  the  way  they  are 
presented  are  in  order.  The  figures  included  here  are  of  two  types.  In  the  first 
type,  we  present  the  radial  cross-sections  of  the  wakes.  Four  cross-sections 
per  wake  spiral.  90  apart  from  each  other,  are  shown:  i.e..  the  topmost  line 
corresponds  to  the  blade  trailing  edge,  the  second  line  is  the  cross-section  of 
the  wake  90c  behind  the  blade,  the  third  line  is  I80c  behind  and  the  fourth  is 
270°  behind  the  trailing  edge.  The  fifth  line  is  360-  or  one  revolution  behind 
the  trailing  edge  and  is  therefore  approximately  on  the  same  plane  as  the 
trailing  edge  (as  are  the  ninth  line,  the  thirteenth  line,  and  so  on).  In  the 
second  class  of  figures,  we  present  the  sectional  lift  coefficient  as  a  function 
of  the  dimensionless  radius  r/R  for  different  time-steps  or  different  cases. 


5.2  Parametric  Study 


A  few  parameters  in  the  formulation  are  empirical.  These  include  on  the 
one  hand,  both  the  initial  size  of  the  vortex-filament  core  and  its  growth  in 
time  (see  Section  3.9).  and  on  the  other  hand,  the  intermediate  wake  model 
(which  depends  on  Cr)  and  the  far  wake  model  which  depends  on  the  estimate 
of  the  flux  through  the  rotor.  The  effects  of  these  parameters  were  studied 
(results  presented  below)  in  order  to  gain  a  better  understanding  of  their 
significance.  In  all  the  results  presented  here,  the  value  of  the  radius  of  the 
vortex  filaments  at  the  trailing  edge  is  set  equal  to  0.3249.  which  corresponds 
to  30ct  of  the  chord.  In  view  of  the  fact  that  the  wake  is  represented  by 
discrete  vortices,  we  estimate  that  this  corresponds  to  60,  c  =  0.2  (the  estimate 
is  based  upon  the  velocity  induced  by  a  finite-thickness  vortex  sheet  and  its 
approximation  with  finite-core  vortex  filaments)  and  seems  to  be  a  barely 
acceptable  upper  estimate  for  the  boundary  layer  thickness.  Lower  sizes  for 
the  vortex  filament  core  appear  desirable;  different  values  were  explored  but 
it  appears  that  this  value  is  necessary  for  the  robustness  of  the  solution. 

An  appreciation  for  the  effects  of  the  size  of  the  vortex-filament  core 
is  obtained  from  the  analysis  of  the  rate  of  growth  of  the  core  size.  Two 
cases  were  considered  for  the  analysis  of  Section  3.9.  In  one,  the  growth  was 
assumed  to  be  such  that  the  core  size  doubled  during  the  first  wake  spiral. 
For  the  thin-wake  analysis  presented  in  Section  3.9  one  obtains  ReE  =  11,845 
(see  Eq.  3.76  with  6i/60  =  2  and  R/c  -  17.5/1.085).  The  other  case  is  that  for 
which  the  core  size  quadruples  in  the  first  blade  revolution:  this  corresponds 
to  Rce  =  2.435.  We  feel  that  these  values  are  realistic  and  bracket  the  actual 
values  for  ReE. 

As  far  as  the  intermediate  and  far-wake  modeling  are  concerned,  the 
following  were  used.  For  all  the  runs,  five  wake  spirals  were  used.  Unless 
stated  otherwise,  only  the  first  wake  spiral  was  treated  as  being  free  whereas 
the  other  four  were  treated  as  the  intermediate  wake  (i.e.,  with  prescribed 
velocity).  This  means  that  at  each  time  step  the  locations  were  recalculated 
for  all  points  of  the  wake,  and  the  contributions  of  all  the  spirals  were  taken 
into  account  in  evaluating  the  velocity  at  the  nodes  of  the  first  spiral.  Of 
course,  all  the  five  spirals  were  taken  into  account  in  solving  the  integral 
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equation. 

The  intensity  of  the  sink  disk  was  evaluated  from  the  velocity  at  the  rotor 
plane  obtained  from  the  actuator  disk  theory  as 

uo  =  fi  R  x C  j  2 

This  value  corresponding  to  an  infinite  blade  analysis  might  be  too  low:  in 
view  of  the  fact  that  the  numerical  results  indicate  a  very  smooth  transition 
between  the  free  portion  and  the  intermediate  portion  of  the  wake,  the  inter¬ 
mediate  wake  expression  might  be  a  better  way  to  estimate  the  flux  through 
the  rotor.  Assessing  u0  from  the  expressions  used  for  the  intermediate  wake 
(see  Landgrebe,  Ref.  12).  we  obtained 

uo  =  1.\5Q>£IR\JCt  i 2 

The  intensity  of  the  disk  is  31  ft/s  with  the  former  model  and  66.7  ft/s  the 
latter.  In  both  cases,  the  disk  was  located  at  2  =  -17.5  ft  (2  =  0  corresponds 
to  the  plane  of  the  rotor)  and  has  a  radius  of  14.0  ft  (this  is  the  estimated 
location  of  the  tip  vortex  after  five  revolutions). 

In  conclusion,  initially,  four  cases  were  analyzed  as  indicated  below  (other 
cases  are  discussed  later):  unless  otherwise  specified,  all  the  cases  are  for 


rotational-flow  formulation. 

Case  #  Strength  6i  /<50 

1  31.0  2 

2  31.0  4 

3  66.7  2 

4  66.7  4 


Case  1  shown  in  Figure  2  is  one  of  the  most  interesting  cases  we  have 
encounter^  and  we  would  classify  it  as  a  ‘borderline’  case  between  stable 
and  unstable  runs.  For  this  reason,  we  present  an  extensive  representation 
of  the  time  history.  The  evolution  of  the  wake  is  shown  by  means  of  the 
plots  of  the  geometry  at  every  six  iterations  (or  time  steps  starting  from 
the  time  step  in  which  the  wake  is  ‘unfrozen.’  see  Section  5.1).  [It  may  be 
pointed  out  that  for  all  the  runs  in  this  report,  the  number  of  time  steps  per 
revolution  is  equal  to  the  number  of  elements  per  spiral,  i.e..  12.  Therefore  six 
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time  steps  correspond  to  one  half  revolution  of  the  rotor. j  One  may  observe 
that  after  two  revolutions,  the  wake  seems  to  go  into  a  limit  cycle,  i.e..  it 
oscillates  radially  w'ith  a  periodic  motion  without  apparently  "going  wild.’ 
The  period  of  the  oscillation  appears  to  be  approximately  two  revolutions 
(compare  iterations  #24  with  #48.  #30  with  #54  and  #36  with  #60).  We 
can  also  describe  the  phenomenon  as  traveling  waves  with  a  wavelength  of 
twice  the  pitch.  It  is  not  clear  to  us  how  these  results  should  be  interpreted 
and  how  close  to  physical  reality  they  are. 

Case  2  (same  disk  strength  but  lower  eddy-Reynolds  number,  i.e.,  6X,  60  = 
4)  is  shown  in  Figure  3  and  indicates  a  much  more  stable  wake.  Only  four 
time  steps  during  the  last  (approximately  the  fourth)  revolution  that  was 
obtained  are  presented.  It  is  apparent  that  the  top  portion  is  completely 
stabilized  and  the  initial  irregularities  are  being  ‘washed  out,’  i.e.,  convected 
by  the  prescribed  velocity  used  for  the  intermediate  wake.  [It  is  apparent 
that  the  model  used  for  the  intermediate  wake  velocity  (from  the  Landgrebe 
model)  is  not  realistic  for  the  tip-vortex  region:  the  tip  vortex  has  a  lower 
velocity  and  therefore  the  vortex  sheet  "catches  up’  with  it,  yielding  an  unre¬ 
alistic  crossing.  However,  the  crossing  does  not  cause  any  problem  since  the 
velocity  of  these  points  are  prescribed.)  It  may  be  worth  noting  that  one  of 
the  points  of  the  wake  spiral  underneath  the  blade  (fifth  line  from  the  top) 
actually  touches  the  blade  itself.  However  this  does  not  seem  to  introduce 
anomalies  in  the  motion  of  the  wake  sheet.  Again,  this  is  unrealistic  and  is 
probably  due  to  the  effect  of  viscosity  (it  may  be  worth  emphasizing  again 
that  in  the  rotational  flow,  we  do  not  consider  the  viscosity  as  an  ‘artificial’ 
one  but  as  an  ‘eddy  viscosity,’  i.e.,  as  a  crude  modeling  of  turbulent  flow). 
On  the  other  hand,  the  sectional  lift  coefficient  presented  in  Figure  4  shows 
an  anomaly  where  the  vortex  touches  the  blade.  In  addition,  it  does  not 
indicate  a  trend  to  form  a  spike  near  the  tip  of  the  blade  as  expected  to 
occur  a  the  wake  rolls  up. 

The  results  for  Case  3  shown  in  Figure  5  are  extremely  stable:  as  men¬ 
tioned  above,  this  corresponds  to  a  disk  strength  of  66.7  ft/s  (obtained  from 
the  flux  through  the  disk  as  evaluated  from  Landgrebe’s  model)  and  a  6,  60 
of  2  (ReE  =  11.845).  It  is  apparent  that  the  increase  in  the  strength  of  the 
sink  has  had  the  effect  of  stabilizing  the  results  obtained  in  Case  1.  On  the 


other  hand,  the  results  appear  to  be  even  smoother  than  in  Case  2  for  which 
we  used  higher  viscosity  (but  lo.ver  sink  strength).  It  may  be  worth  noting 
that  increasing  the  strength  to  100  ft/s  (Case  3A.  not  shown  here)  does  not 
seem  to  affect  the  wake  motion  within  plotting  accuracy.  The  sectional  lift 
distribution  is  shown  in  Figure  6  and  again  no  spike  seems  to  be  present 
indicating  that  the  solution  obtained  for  the  wake  geometry  might  be  too 
smooth. 

In  Case  3B  which  is  the  same  as  Case  3  except  that  the  potential-flow 
formulation  is  used  instead  of  the  rotational-flow  formulation  (as  pointed 
out  above,  in  this  case  the  viscosity  may  be  introduced  only  as  artificial 
viscosity).  The  results,  shown  in  Figure  7.  are  extremely  similar  to  those 
of  Figure  5.  Similar  analogy  was  encountered  in  all  the  runs  that  we  made 
for  the  potential-flow  formulation,  indicating  the  computational  equivalence 
(not  the  conceptual  one)  of  the  two  methods. 

Next,  consider  Case  4.  shown  in  Figure  8.  The  effect  of  the  increase 
in  viscosity  over  Case  3  appears  to  be  a  sharper  turning  of  the  wake  vortex 
sheet.  In  fact,  this  trend  is  very  similar  to  that  of  Case  2  (compare  Figure  5D 
to  Figure  3C,  both  corresponding  to  the  48th  time  step).  The  effect  of  the 
higher  suction  in  Case  4  has  the  predictable  effect  of  increasing  the  distance 
of  the  tip  vortex  corresponding  to  the  first  spiral,  from  the  blade  (in  Case 
2,  the  tip  vortex  of  the  first  spiral,  fifth  line  from  the  top.  was  touching  the 
blade).  Surprisingly  enough.  Case  4  has  a  higher  value  of  the  sectional  or 
local  lift  coefficient  (see  Fig.  9)  than  Case  2  (Figure  4),  where  the  vortex 
was  touching  the  blade. 

In  all  of  the  above  results,  only  the  first  spiral  is  free.  In  order  to  assess 
the  effect  of  using  only  one  free  wake,  consider  Case  5.  which  is  similar  to 
Case  4.  but  with  the  first  two  wake  spirals  being  free  (with  velocity  calculated 
at  each  time  step),  while  the  intermediate  wake  is  composed  of  the  last  three 
wake  spirals.  This  is  the  longest  run  that  we  made  eighty  time  steps, 
or  almost  seven  revolutions.  The  transient  behaviour  of  the  wake,  shown 
in  Figure  10.  is  similar  to  that  presented  in  Ref.  1.  The  wake  geometry 
becomes  very  messy  in  the  region  of  the  second  spiral  around  the  24th  time 
step.  Around  the  36th  time  step,  the  first  two  spirals  are  smoother  and 
eventually  the  solution  converges  (which  was  not  the  case  in  Ref.  1  because 
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of  the  lack  of  intermediate  and  far-wake  models).  The  time  history  of  the 
loading  distribution  is  shown  in  Figure  11.  The  spike  is  more  pronounced 
than  in  Case  4  (one  free  spiral). 

As  apparent  from  Figure  12.  comparison  with  the  results  of  Rao  and 
Schatzle  (Ref.  7)  and  Morino.  Kaprielian  and  Sipcic  (Ref.  1)  show  that 
Case  4  (one  free  spiral)  and  Case  5  (tw’o  free  spirals)  are  in  good  agreement 
with  the  results  of  Rao  and  Schatzle  (who  show  excellent  agreement  with 
experimental  results  for  the  case  of  a  four-bladed  rotor,  see  Section  5.3).  It 
should  be  noted  that  Rao  and  Schatzle  obtained  results  by  a  generalized 
wake  analysis  (using  Landgrebe  model  for  the  wake  geometry)  whereas  our 
results  are  obtained  using  a  free-wake  analysis.  Also  shown  in  the  figure 
are  the  results  of  our  own  generalized  wake  analysis  (also  using  a  modified 
Landgrebe  wake)  which  are  in  good  agreement  with  the  two-free-spiral  free- 
wake  analysis. 

5.3  Multi-Bladed  Rotor 

The  results  presented  in  Section  5.2  are  all  for  a  single  bladed  rotor.  In 
order  to  have  an  assessment  on  more  realistic  configurations,  two-bladed  and 
four-bladed  rotors  have  been  investigated.  Results  from  a  prescribed  wake 
analysis  are  given  for  a  two-bladed  and  a  four-bladed  rotor,  while  the  free 
wake  analysis  is  limited  to  the  case  of  a  two-blade  case. 

Consider  Case  6.  i.e.,  the  free-wake  analysis  for  a  two-bladed  rotor.  The 
data  used  for  Case  6  are  very  similar  to  that  of  Case  4  with  the  number 
of  blades  increased  to  two:  the  number  of  free  spirals  is  retained  as  one. 
The  wake  geometry  after  45  iterations  for  Case  6  are  shown  in  Figure  13. 
One  may  observe  the  strong  similarity  between  this  and  the  wake  of  Case 
4  for  the  corresponding  time-step.  If  the  computations  are  continued  for 
additional  time-steps,  we  expect  the  wake  to  stabilize  as  for  the  single-bladed 
rotor.  The  sectional  lift  distribution  is  shown  in  Figure  14  along  with  results 
from  our  own  generalized  wake  analysis.  Again,  the  results  are  in  good 
agreement.  This  gave  us  sufficient  confidence  in  our  own  generalized-wake 
analysis  to  use  it  for  an  indirect  assessment  of  the  method  (the  free-wake 
analysis  for  the  four-bladed  rotor  has  not  yet  been  attempted):  the  results 
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for  the  generalized-wake  analysis  for  a  four-bladed  rotor  are  compared  in 
Figure  14A  with  the  numerical  results  of  Rao  and  Schatzle  (Ref.  7)  and 
the  flight  test  data  presented  by  Bartsch  (Ref.  16)  for  a  XH-51A  helicopter 
rotor. 

5.4  Ground  Effect 

The  final  problem  considered  here  is  the  case  of  a  single-bladed  rotor  in 
ground  effect.  Contrary  to  our  expectations,  as  we  lowered  the  rotor-plane 
to  the  ground,  the  problems  we  encountered  were  not  in  the  region  near 
the  ground,  but  were  in  the  region  near  the  rotor.  We  were  able  to  explain 
this  with  the  fact  that  as  the  image  became  closer  to  the  rotor,  it  had  the 
effect  of  a  source  (negative  sink):  the  spirals  got  closer  and  the  effect  of  the 
concentrated  vortex  discretization  became  more  pronounced.  [As  a  point 
gets  close  to  a  vortex  sheet,  the  velocity  induced  is  tangential  to  the  sheet. 
On  the  other  hand,  if  a  point  gets  close  to  a  series  of  of  concentrated  vortices 
simulating  the  vortex  sheet,  the  velocity  is  tangential  to  the  sheet  only  if 
the  point  approaches  a  vortex  along  the  normal  to  the  vortex  sheet  through 
that  vortex.  Very  strong  normal  components  arise  from  any  deviation  from 
that  line.]  Therefore,  we  limited  our  analysis  to  a  rotor  which  is  sufficiently 
close  to  the  ground  to  experience  ground  effect,  but  not  so  close  that  the 
above  described  effect  of  concentrated  vortex  discretization  would  require 
extremely  high  viscosity  to  'stabilize-  the  solution. 

The  results  presented  here.  Case  7,  are  for  the  same  problem  as  Case  4, 
except  for  ground  effect.  The  initial  wake  geometry  is  the  same  as  that  for  the 
out-of  ground-effect  test  cases.  The  prescribed  velocity  in  the  intermediate 
wake  is  such  that  the  wake  generated  would  have  to  go  to  zero  exponentially 
with  v,  the  azimuth  angle,  whereas  the  radius  would  go  exponentially  to 
infinity  (see  Ref.  38  for  details).  The  far-wake  model  is  a  distribution  of 
sinks  over  a  cylinder  of  20  ft  radius.  The  intensity  is  obtained  for  the  flow 
through  the  rotor  as  evaluated  from  the  actuator  disk  theory. 

The  wake  geometry  at  different  time-steps  is  showrn  in  Figure  15.  The 
results  are  extremely  smooth:  this  appears  to  be  more  the  exception  than 
the  rule  for  the  ground-effect  problem.  Figure  16  shows  the  sectional  lift  dis- 


tribution.  These  results  are  presented  here  more  in  terms  of  addressing  the 
issue  of  the  stability  of  the  wake  than  for  an  actual  comparison  with  existing 
results.  It  is  acknowledged  that  considerable  additional  work  is  required. 
However,  we  feel  that  meaningful  comparisons  with  exisiting  data  can  be 
made  only  after  the  numerical  phenomenon  described  above  (velocities  nor¬ 
mal  to  the  vortex  sheet  introduced  by  the  concentrated  vortex  discretization) 
is  eliminated  by  introducing  a  first-order  formulation. 


SECTION  6 


CONCLUSIONS 


Two  integral  equation  methods  for  the  free-wake  analysis  of  helicopter 
rotors  (for  potential  and  rotational  flows,  respectively)  have  been  presented. 
The  rotational-flow  formulation  is  based  upon  Helmholtz  scalar/ vector-po¬ 
tential  decomposition.  The  advantages  of  the  rotational-flow  formulation 
over  the  potential-flow  formulation  have  been  discussed.  The  numerical 
equivalence  of  the  two  methods  has  been  demonstrated.  It  should  be  noted 
again  that,  whereas  in  the  potential-flow  formulation  viscosity  may  be  in¬ 
troduced  only  as  artificial  viscosity,  in  the  rotational-flow  formulation  the 
presence  of  viscosity  (or  eddy-viscosity)  is  consistent  with  the  correspond¬ 
ing  assumptions.  Therefore,  this  formulation  is  applicable  to  the  solurion  of 
time-averaged  Navier-Stokes  equations:  in  particular,  a  simplified  thin-wake 
analysis  with  an  elementary  eddy-viscositv  model  for  turbulence  was  used 
for  the  numerical  applications. 

Results  obtained  indicate  that,  for  appropriate  values  of  empirical  param¬ 
eters  (e.g.,  thickness  of  the  wake  at  the  trailing  edge,  eddy- Reynolds  number, 
intensity  of  the  far  wake  sink  disk),  the  solution  reaches  steady-state,  with 
the  wake  converging  to  a  smooth  geometry  and  the  sectional  lift  distribution 
in  good  agreement  with  the  numerical  results  of  Rao  and  Schatzle  [Ref.  7] 
and.  indirectly,  with  the  experimental  ones  of  Bartsch  [Ref.  16] .  However, 
a  sensitivity  analysis  over  the  effects  of  these  paramenters  indicates  that  the 
results  were  not  always  as  good.  For  some  values  of  these  parameters,  the 
solution  may  reach  a  steady  state  but  the  sectional  lift  distribution  does  not 
necessarily  show  a  marked  improvement  as  the  wake  rolls  up.  For  other  val¬ 
ues  of  these  parameters,  the  solution  may  be  completely  unstable.  The  main 
reason  for  this  sensitivity  to  the  values  of  the  empirical  paramenter  is  at¬ 
tributed  to  the  fact  that  concentrated  vortices  are  used  for  the  discretization 


of  the  formulation  (see  Section  5.4).  Ground  effect  results  are  also  included, 
and  show  that  stable  results  with  the  proper  trend  may  be  obtained. 

Most  of  the  computations  were  performed  on  a  VAX-780  computer.  Case 
7,  for  ground  effect,  took  10  hours  and  58  minutes  to  perform  72  iterations 
(in  addition  to  62  time  steps  with  frozen  wake).  It  takes  one  hour  to  generate 
12  time  steps  for  a  single-bladed  rotor.  Some  runs  were  made  on  an  IBM 
3081.  The  longest  run  (Case  5,  for  single-bladed  rotor  with  two  free  spirals) 
took  2  hours  and  2  minutes  to  perform  80  iterations.  It  took  one  hour  and 
18  minutes  for  Case  4,  for  single-bladed  rotor  with  one  free  spiral,  also  for  80 
iterations.  The  potential-flow  formulation  requires  approximately  the  same 
amount  of  computer  time  as  the  rotational-flow  formulation  (a  significantly 
shorter  time  than  that  encountered  by  Morino.  Kaprielian  and  Sipcic,  Ref. 
1).  In  summary,  we  may  conclude  that  the  new  formulation  presented  here 
(for  rotational  flows)  has  much  broader  applicability  than  the  formulation  for 
potential  flows,  while  requiring  approximately  the  same  amount  of  computer 
time. 

In  conclusion,  we  believe  that  the  investigation  was  completed  in  terms  of 
the  original  goals  of  the  project,  i.e.,  development  of  a  rotational-flow  formu¬ 
lation  as  efficient  as  the  potential-flow  formulation,  but  free  of  the  inherent 
limitations  of  the  latter.  However,  the  investigation  should  not  be  consid¬ 
ered  as  completed,  since  it  has  not  addressed  in  a  fully  satisfactory  fashion 
all  the  issues  encountered.  In  particular,  in  terms  of  recommendations  for 
future  work,  we  believe  that  the  most  important  item  to  be  addressed  is  the 
first-order  formulation  for  the  wake  (as  outlined  in  Section  4).  As  discussed 
in  Section  5.4,  the  main  problem  encountered  in  this  effort  has  to  do  with  the 
presence  of  induced  velocities  normal  to  the  wake  sheet:  this  is  attributable 
exclusively  to  the  zeroth-order  formulation.  It  may  be  noted  that  a  first- 
order  formulation  (as  outlined  in  Section  4)  has  been  implemented,  under  a 
separate  effort,  for  the  problem  of  an  airplane  wing  [Ref.  39].  The  incorpo¬ 
ration  of  this  formulation  into  the  code  for  free-wake  analysis  of  helicopter 
rotors  is  not  expected  to  present  any  significant  difficulty,  but  the  issue  of 
numerical  instabilities  needs  to  be  addressed  again. 

A  second  item  which  requires  further  investigation  is  the  modelling  for 
the  intermediate  and  far  wakes.  We  do  not  believe  that  we  have  addressed 


this  model  are  now  given  as  input  data:  self-adjusting  models  should  be 
incorporated  into  the  formulation.  Similar  comments  should  be  made  about 
the  turbulence  modelling  and  the  boundary  layer  analysis:  while  w’e  believe 
that  we  have  shown  that  actual  turbulence  modelling  and  boundary  layer 
analysis  can  be  incorporated  in  the  formulation,  a  considerable  amount  of 
work  is  needed  to  accomplish  that,  even  for  the  relatively  simple  problem  of 
attached  flow’s. 

After  these  items  have  been  addressed,  more  extensive  convergence  anal¬ 
ysis  and  comparison  of  the  numerical  results  with  existing  experimental  re¬ 
sults  should  be  undertaken,  especially  as  far  as  the  ground  effect  problem  is 
concerned. 

Finally,  it  should  be  emphasized  that  all  of  the  results  were  obtained 
through  a  time-accurate,  step-by-step  solution.  Therefore,  the  formulation 
is  applicable,  at  least  in  theory,  to  the  solution  of  more  general  time-domain 
problems.  The  accuracy,  robustness  and  efficiency  of  the  algorithm  should  be 
addressed  by  applying  it,  for  instance,  to  helicopter  rotors  in  forward  flight. 
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Figure  2D.  Wake  Geometry  for  Case  1 
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Figure  2E.  WaXe  Geometry  for  Case  1 
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Figure  2K.  Wake  Geometry  for  Case  1 
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Figure  3C.  Wake  Geometry  for  Case  2 
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Figure  3D.  Wake  Geometry  for  Case  2 
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Figure  6.  Sectional  Lift  Distribution  for  Case 
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Figure  7B.  Wake  Geometry  for  Case  3a 
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Figure  8a.  Wake  Geometry  for  Case  4 
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Figure  8D.  Wake  Geometry  for  Case  4 
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Figure  8E.  Wake  Geometry  for  Case  A 
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Figure  10G.  Wake  Geometry  for  Case  5 
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Figure  10H.  Wake  Geometry  for  Case  5 
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Figure  101.  Wake  Geometry  for  Case  5 
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SINGLE  BLADED  ROTOR 


Figure  12.  Comparison  of  present  results  for  one-bladed  rotor  (free  wake 
analysis  of  cases  A  and  5  and  generalized  wake  analysis)  with  those  of  Rao 
and  Schatzle  (Ref.  7)  and  Morino,  Kaprielian,  and  Sipic  (Ref.l). 
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13D.  Wake  Geometry  for  Case  6 
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Figure  13E.  Wake  Geometry  for  Case 
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Figure  13J.  Wake  Geometry  for  Case  6 


Figure  13K.  Wake  Geometry  for  Case  6 
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Figure  13L.  Wake  Geometry  for  Case  6 
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Figure  14.  Current  results  for  two-bladed  rotor  (free-wake  and 
generalized-wake  analysis). 


FOUR  BLADED  ROTOR 


Figure  14A. Conpari son  of  present  results  for  four-biaded  rotor  (generalized 
wake  analysis)  with  those  of  Rao  and  Schatzle  (Ref.  7)  and  the  experiaents 
of  Bartsch  (Ref.  16). 
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Figure  15B.  Wake  Geometry  for  Case  7 
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Figure  15D.  Wake  Geometry  for  Case  7 
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Figure  15K.  Wake  Geometry  for  Case  7 
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Figure  16.  Sectional  Lift  Distribution  for  Case 


